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1.  INTRODUCTION 


In  Ref.  [  1  ]  new  methods  for  structural  weight  minimization 
were  presented,  which  were  restricted  to  finite  element  models 
made  up  of  bars  and  membranes.  The  only  constraints  imposed  in 
the  optimization  problem  were  concerned  with  static  stress  and 
displacement  limitations.  The  main  goal  of  the  present  work 
is  to  extend  the  mixed  method  and  the  generalized  optimality 
criterion  set  forth  in  Ref,  [  1  ]  to  deal  with  more  complex  pro¬ 
blems  involving  frequency  and  buckling  constraints  as  well  as 
bending  elements. 

It  has  been  shown  in  Ref.  [  1  ]  that  optimality  criteria  and 
mathematical  programming  approaches  to  structural  optimization 
have  now  reached  a  stage  where  they  use  the  same  basic  concepts. 
Both  approaches  convert  the  primary  minimization  problem  into 
a  sequence  of  explicit  approximate  problems.  To  this  end  con¬ 
ventional  optimality  criteria  techniques  are  based  upon  physical 

understanding  of  the  phenomena  involved  when  redesign  takes  place 
* 

(load  pathes,  stress  redistribution,  etc...).  Another  intuitive 
guideline  is  employed  for  deriving  recursive  redesign  formulas. 
However  these  optimality  criteria  techniques  can  be  justified 
from  a  mathematical  programming  point  of  view,  yielding  therefore 
a  unified  approach  to  structural  optimization.  This  unified  ap¬ 
proach  will  be  reviewed  in  chapter  2  with  emphasis  on  the  primal 
or  dual  character  of  the  solution  scheme  used  to  solve  each  ex¬ 
plicit  prob 1 em . 


As  in  Ref.  (  1  )  the  design  variables  will  be  taken  as  the 
transverse  sizes  of  the  s  iuctural  members,  namely,  the  cross- 
sectional  areas  of  bar  at  '  beam  elements  and  the  thicknesses  of 
membrane,  plate  and  flat  shell  elements.  The  objective  func¬ 
tion  to  be  minimized  is  defined  as  a  linear  function  of  the  mem¬ 
ber  sizes  a  ^  : 


W 


n 


i*i 


a  . 
i 


(  I  .  I  ) 


where  is  a  geometrical  factor  such  that  the  product  a^  is 
the  volume  of  the  it'1  element  (length  of  aone-dimensional  element; 
area  of  a  two-dimensional  element).  P.  denotes  a  scalar  quantity 
associated  with  the  i ^  element.  In  most  cases,  it  is  simply  the 
specific  weight,  in  which  case  the  objective  function  posed  by 
Eq .  (1.1)  represents  the  weight  of  the  structure. 

■^The  structural  optimization  problem  considered  in  this  work 
consists  in  minimizing  the  objective  function-,  embodied  in  Eq  . 

(  I  .  1 )  -?sub  j  ect  to  a  wide  variety  of  constraints.  The  behavior 
constraints  impose  limitations  on  the  quantities  describing  the 
structural  response  :  static  stresses  and  displacements  under  mul¬ 
tiple  loading  conditions,  natural  frequencies  and  linear  buckling 
load  factors.  In  addition,  side  constraints  are  directly  pres¬ 
cribed  on  the  design  variables  themselves., 

1 . 1  Side  constraints 

For  various  reasons  it  is  generally  necessary  to  impose  lower 
and  upper  limits  on  the  design  variables,  which  constitutes  the 
so-called  side  constraints  : 

a .  <  a .  <  a .  (1.2) 

—ill 

Note  that  throughout  this  work,  a  lower  bar  denotes  a  minimum 
allowable  limit,  while  an  upper  bar  represents  a  maximum  allowa¬ 
ble  limit.  Although  the  side  constraints  embodied  in  Eq .  (1.2) 

do  not  correspond  really  to  failure  modes  of  the  structure,  they 
are  nevertheless  important  in  the  problem  statement,  because  they 
reflect  fabrication  and  analysis  validity  considerations. 

It  is  also  desirable  to  impose  several  structural  members  to 
be  described  by  the  same  design  variable.  Therefore  the  finite 
elements  are  subdivided  into  some  preselected  groups,  so  that  one 
independent  design  variable  is  sufficient  to  control  the  size  of 
all  finite  elements  pertaining  to  the  same  linking  group.  Design 
variable  linking  can  be  introduced  in  the  form  of  equality  cons- 
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traints  on  the  member  sizes  : 

a  =  G  a '  (1.3) 

where  a  (dimension  n)  is  the  full  vector  of  the  member  sizes, 
a*  (dimension  n')  is  the  reduced  vector  of  independent  design 
variables  and  G  (dimension  n  x  n')  is  a  boolean  matrix  (0  or  1 
elements).  Design  variable  linking  reduces  the  number  of  design 
variables  (n'  <  n)  and  it  facilitates  imposition  of  constraints 
based  on  symmetry,  fabrication  requirements  and  cost  considera¬ 
tions  regarding  the  number  of  parts  to  be  assembled.  In  chapter 
2  through  6  no  attention  will  be  paid  to  the  treatment  of  the 
linking  constraints  expressed  in  Eq.  (1.3),  however  chapter  7 
will  indicate  how  to  handle  them  in  a  simple  way. 

1.2.  Stress  and  displacement  constraints 

A  significant  class  of  structural  optimization  problems  is 
concerned  with  constraints  on  the  static  structural  response  under 
a  set  of  loading  conditions.  Strength  constraints  are  usually 
taken  into  account  by  p re s cr ib ing a  1 lowab le  limits  on  the  stresses 
in  the  various  members  : 

°kJl  <  °k  k  -  I  ,n  (1.4) 

£  =  1  ,c 

where  o^  denotes  an  adequate  representation  of  the  stress  state 
in  the  kc^  element  under  the  £t^1  loading  condition  and  o^  is  the 
allowable  stress  limit  in  this  particular  element.  For  stress 
members  independent  tension  and  compression  allowables  can  be 
specified.  In  shear  panel  and  isotropic  membrane  elements,  where 
multiaxial  stress  states  exist,  strength  constraints  are  introduced 
by  limiting  the  value  of  an  equivalent  stress,  using  for  example 
the  well  known  von  Mises  criterion.  In  the  orthotropic  membrane 
elements  employed  to  model  fiber  composite  lamina^ several  strength 
failure  criteria  can  be  adopted,  for  example,  the  Tsai-Hill  cri¬ 
terion. 
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The  displacement  constraints  considered  in  Ref.  [  1  ]  are  » 
defined  as  upper  bounds  on  some  linear  combinations  of  the  dis¬ 
placement  degrees  of  freedom  q  used  in  the  finite  element  model 


.T 

bj  q 


(1.5) 


where  bj  is  a  vector  of  constants.  The  linear  combination  of 
displacements  u.  is  called  a  flexibility.  The  displacement  or 
flexibility  constraints  read  thus  as  follows  : 


u  .  <  u  . 

J  SL  J 


(  1  .6) 


where  Uj^  denotes  the  value  of  Uj  under  the  l  loading  condition 
and  Uj  is  the  upper  limit.  Note  that  this  treatment  of  displace¬ 
ment  constraints  includes  the  usual  nodal  displacement  constraints, 
the  relative  displacement  constraints,  the  slope  constraints,  etc... 

In  the  displacement  finite  element  method  employed  in  Ref.  I  1  ]  , 
as  well  as  in  the  present  work,  the  structural  analysis  needed  for 
evaluating  the  stress  and  flexibility  constraints  embodied  in  Eqs. 
(1.4,  1.6)  is  performed  by  solving  the  systems  of  linear  equations: 


(1.7) 


where  K  is  the  structural  stiffness  matrix  and  g  represents  the 
forces  acting  in  the  i  loading  condition.  Once  the  displacement 
vector  q  is  known,  the  stresses  in  the  various  members  are  eva¬ 


luated  from  the  equations  : 


Tk  q* 


(  1.8) 


where  Tk  is  the  stress  matrix  for  the  k  element.  Note  that  in 
a  membrane  element  the  stress  matrix  is  made  up  of  three  rows, 
while  in  a  bar  element,  it  reduces  to  one  row.  The  class  of  fi¬ 
nite  element  models  studied  in  Ref.  [  1  ]  was  restricted  to  thin- 
walled  structures  (assembing  of  bars  and  membranes)  for  which  the 


stress  matrices  are  independent  of  the  design  variables  and  the 
stiffness  matrix  exhibits  a  linear  form  in  the  design  variables  : 

n 

K  =  . E .  a .  K  .  (1.9) 

1=1  1  l 

where  represents  the  stiffness  matrix  of  the  ith  element  when 
a .  =  I  . 

l 

Chapter  3  will  be  devoted  to  recalling  the  main  results  ob¬ 
tained  in  Ref.  [I  ] .  An  important  result  is  that  there  exists 
two  distinct  but  equivalent  procedures  for  generating  first  order 
explicit  approximations  of  the  stress  and  displacement  constraints. 
Both  procedures  require  that  a  certain  number  of  additional  loa¬ 
ding  cases  be  treated  in  the  structural  analysis  phase.  The  gene¬ 
ralized  optimality  criterion  set  forth  in  Ref.  [  I  ]  employs  first 
order  approximation  for  all  the  stress  and  displacement  constraints. 
It  can  be  written  in  terms  of  virtual  strain  energy  densities  in 
the  structural  members.  In  conventional  optimality  criteria  tech¬ 
niques,  only  the  displacement  constraints  are  replaced  with  first 
order  approximations,  while  the  stress  constraints  are  transformed 
into  minimum  size  limits,  which  corresponds  to  adopting  zero  order 
explicit  approximations  (fully  stressed  design  concept).  There¬ 
fore  hybrid  optimality  criteria  can  be  defined  as  those  based  on 
mixed  order  approximations  of  the  stress  constraints.  Note  that 
the  results  stated  in  Ref.  [  1  ]  were  concerned  with  applied  loads 
that  are  independent  of  the  design  variables.  They  will  be  ex¬ 
tended  in  this  work  to  the  case  where  the  load  vectors  g^  appearing 
in  Eq.  (1.7)  depend  linearly  on  the  design  variables  (thermal  ef¬ 
fects,  body  forces,  etc...). 

1.3.  Frequency  constraints 

The  natural  frequencies  of  a  structure  are  often  constrained 
to  reside  within  some  given  critical  intervals  in  order  to  insure 
adequate  control  of  integrated  mechanical  systems,  or  to  prevent 
resonance  phenomena,  etc...  The  frequency  constraints  will  there¬ 
fore  be  expressed  as  follows  : 
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(1)  .  <  0)  .  <  0)  . 
-J  J  J 


j  =  1  ,m 


( I . 10) 


where  the  “ j ' s  represent  the  eigenf requencies  ordered  in  an  in¬ 
creasing  sequence.  Most  often  only  the  fundamental  natural 
frequency  Wj  is  required  to  be  larger  than  an  allowable  limit, 
in  which  case  Eq .  (I. JO)  reduces  to  a  single  behavior  constraint. 

Note  however  that  when  the  optimization  process  progresses,  a 
change  might  occur  in  the  lowest  vibration  mode,  so  that  it  is 
more  reliable  to  always  impose  multiple  frequency  constraints, 
at  least  on  the  two  first  eigenmodes. 


Instead  of  solving  systems  of  linear  equations,  structural 
analysis  now  consists  of  solving  an  eigenproblem 


Kq.  -  ut  M  q  .  =  0 
J  J  J 


d.ll) 


where  u>  j  is  the  j frequency,  q^  denotes  the  corresponding  mo¬ 
dal  displacement  vector  and  M  represents  the  structural  mass 
matrix.  It  will  be  seen  in  chapter  4  that  the  frequency  cons¬ 
traints  can  be  treated  exactly  like  the  stress  and  displacement 
constraints,  by  using  first  order  Taylor  series  expansion  with 
respect  to  the  reciprocal  design  variables.  The  resulting  ex¬ 
plicit  problem  can  still  be  solved  by  recoursing  to  primal  or 
dual  algorithms,  which  corresponds  to  applying  the  concept  of 
mixed  method  or  that  of  generalized  optimality  criterion,  respec¬ 
tively.  In  particular  the  optimality  criterion  can  be  stated  in 
terms  of  elastic  and  kinetic  energy  densities  in  the  structural 
membc  r s . 

1 . 4  Buckling  constraints 


In  many  design  projects,  stability  considerations  play  an 
important  role  in  checking  off  the  possible  failure  modes  of 
the  structure.  For  numerical  optimization  purposes  it  is  con¬ 
venient  to  distinguish  two  d i f f eren t ty pe s  of  instability  pheno¬ 
mena  :  local  instability  of  an  individual  structural  component 
and  global  instability  of  the  entire  structure.  Local  instability 
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can  often  be  controlled  by  defining  adequate  "equivalent"  stresses 
in  the  strength  constraints  posed  in  Kq  .  (1.4)  and/or  by  reducing 
the  compressive  allowable  stresses  in  axial  members  and  panels. 

For  example  Euler  buckling  constraint  in  a  tubular  member  with 
specified  mean  radius  r.  Young's  modulus  E  and  length  Si  is  ob¬ 
tained  by  replacing  the  compression  allowable  limit  with 

2 

-  2  F  r 

a  =  ir  (1.12) 

2J r 

On  the  other  hand  global  instability  constraints  depend  in  a 
complex  way  on  all  structural  components.  In  the  context  of  a  linear 
analysis,  global  instability  is  taken  into  consideration  through 
an  eigenproblem  similar  to  that  defining  the  vibration  modes  of 
a  structure  : 


K  -  X.  K0  -  0 


0.13) 


where  K^,  is  the  geometric  stiffness  matrix.  The  e  i  gens  o  1  u  t  i  on  s 
of  this  problem  yield  the  buckling  mode  vectors  q^  and  the  as¬ 
sociated  buckling  load  factors  X^.  When  several  distinct  loading 
conditions  are  applied  to  the  structure,  buckling  constraints 
must  be  prescribed  under  any  of  the  loading  conditions,  which 
implies  solving  several  e igenprob lems  at  each  redesign  stage. 

So  the  buckling  constraints  are  generally  imposed  on  several 
critical  loads  corresponding  to  either  different  buckling  modes 
or  to  different  loading  conditions  : 


X  .  <  X  .  <  X  . 
-J  J  J 


(i  •  i  M 


Strictly  speaking  only  the  lower  bound  X^ .  is  significant  for 
representing  global  buckling  contraints,  with  usually  X^  =  I  .  However  it 

is  sometimes  interesting  to  include  multiple  constraints  with 
upper  bounds  X  ^  to  make  possible  the  separation  of  buckling  modes 
that  are  close  to  each  other.  This  is  useful  for  reducing  the 
structural  sensibility  to  initial  imperfections. 


The  buckling  constraints  stated  in  Eq .  (1.14)  will  be  con- 


sidered  in  chapter  5,  using  again  the  concepts  of  mixed  method 
and  of  generalized  optimality  criterion.  The  optimal  design 
will  be  characterized  in  terms  of  strain  energy  densities  in 
the  structural  members  for  the  critical  buckling  modes. 


The  structural  optimization  problem  studied  in  the  present 
work  consists  thus  in  minimizing  the  objective  function  defined 
in  Eq.  (1.1),  subject  to  the  side  constraints  given  in  Eqs.(l.2, 
1.3)  and  to  the  various  behavior  constraints  embodied  in  Eqs  . 

(1.4,  1.6,  1.10,  1.14).  This  problem  is  examined  initially  by 

restricting  its  formulation  to  structural  models  made  up  of  bar, 
membrane  and  shear  panel  elements,  which  are  quite  adequate  for 
idealizing  thin  walled  structures  subject  mainly  to  extension 
loads  (chapters  2  through  5).  Subsequently  chapter  6  will  be 
concerned  with  structural  models  that  are  capable  of  carrying 
flexural  forces,  such  as  beam,  plate  and  shell  elements.  It  is 
shown  that  high  quality  explicit  approximation  of  the  behavior 
constraints  can  still  be  generated  using  first  order  Taylor  se¬ 
ries  expansion,  provided  that  adequate  intermediate  variables 
are  selected.  The  idea  of  generalized  optimality  criterion 
remains  fully  valid  and  it  keeps  its  interpretation  in  terms  of 
energy  densities  in  the  structural  members. 

In  chapter  7  some  indications  will  be  given  about  the  com¬ 
puter  implementation  of  the  optimization  strategies  proposed  in 
Ref.  [  1  1  and  continued  in  the  present  work.  Finally  some  illus¬ 
trative  examples  will  be  presented  in  chapter  8.  Based  on  the 
numerical  results  obtained,  it  can  be  concluded  that  the  concepts 
developed  in  this  work  are  rather  general.  They  lead  to  efficient 
structural  optimization  methods,  which  are  capable  of  dealing  with 
a  wide  variety  of  behavior  constraints  and  of  structural  models. 
The  approach  presented  converts  the  primary  problem  into  a  se¬ 
quence  of  simple  explicit  problems.  The  explicit  problem  state¬ 
ment  is  obtained  by  linearizing  the  behavior  constraints  with 


. . . .  Ml  .  . . .  * _ _ _  .  ......  .. 


respect  to  appropriate  intermediate  variables.  Therefore  the 
method  requires  numerical  evaluation  of  the  constraint  gradients 
at  each  redesign  stage.  Using  a  dual  solution  scheme  for  sol¬ 
ving  each  explicit  problem  leads  to  generalization  of  the  con¬ 
ventional  optimality  criteria  techniques,  with  comparable  ef¬ 
ficiency  but  increased  reliability.  Solving  partially  each 
explicit  problem,  by  employing  a  primal  solution  scheme,  yields 
the  mixed  method,  which  facilitates  control  over  convergence 
of  the  overall  optimization  process. 
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2.  MIXED  METHOD  AND  GENERALIZED  OPTIMALITY  CRITERION 

This  chapter  summarizes  some  theoritical  results  obtained  in  a 
previous  work  [  1,2  ]  .  It  is  shown  that  a  powerful  approach 
to  structural  optimization  has  now  emerged,  which  consists  in 
replacing  the  primary  problem  with  a  sequence  of  simple  expli¬ 
cit  problems.  In  the  optimality  criteria  and  mathematical  pro¬ 
gramming  approaches,  the  behavior  constraints  are  approximated 
using  respectively  virtual  load  considerations  and  linearization 
wi  tli  respect  to  the  reciprocal  design  variables. 

An  attractive  strategy  is  to  solve  partially  each  explicit 
problem,  using  a  primal  solution  scheme,  before  reanalyzing  the 
structure  and  updating  the  approximate  problem  statement.  This 
process  facilitates  generation  of  a  sequence  of  steadily  impro¬ 
ved  feasible  designs.  This  approach  can  be  interpreted  as  a 
mixed  primal-linearization  method  and  it  introduces  an  interes¬ 
ting  possibility  of  controlling  the  convergence  of  the  optimi¬ 
zation  process. 

An  alternative  approach  is  to  recognize  that  the  explicit 
but  approximate  problem  statement  is  of  such  high  quality  that 
it  can  be  solved  exactly,  rather  than  partially.  Such  an  exact 
solution  can  be  efficiently  generated  by  using  a  dual  solution 
scheme.  This  approach  is  closely  related  to  pure  linearization 
methods  in  mathematical  programming,  but  it  can  also  be  viewed 
as  a  rigorous  generalization  of  the  conventional  optimality  cri¬ 
teria  techniques. 


2.1.  Fundamental  concepts 


For  presenting  the  fundamental  concepts  used  in  the  methods 
formulation,  the  optimization  problem  stated  in  chapter  1  is 
written  in  the  mathematical  form  : 

n 

minimize  W  =  .1 ,  P.l.a.  (2.1) 

i=  1  l  l  l 

subject  to 

hj(a)^0  j=!,m 

a  .  >  a  .  >  a  .  i  =  I  ,  n 

li—i 

where  a^  denote  the  n  design  variables  which  correspond  to 
member  sizes  of  either  individual  finite  elements  or,  more 
often  of  some  preselected  groups  of  finite  elements  (see  sec¬ 
tion  1.1  ).  The  objective  function  to  be  minimized  is  the 
structural  weight.  It  is  a  linear  function  of  the  design 
variables  a..  All  the  behavior  constraints  discussed  in  chap¬ 
ter  1  (see  Eq s .  1.6 , 1.14)  are  gathered  in  the  inequalities  ex¬ 
pressed  in  Eq.  (2.2).  They  represent  limitations  on  quantities 
describing  the  structural  response,  for  example,  the  stresses 
and  the  displacements  under  multiple  static  loading  cases,  the 
natural  frequencies,  the  buckling  loads,  etc...  The  design  va¬ 
riables  are  also  subjected  to  side  constraints  that  prevent 
them  from  becoming  too  small  or  too  large  (Eq.  2.3).  Note  that 
the  design  variables  are  assumed  to  vary  continuously  in  Eqs. 

(2. 1-2.3).  However  di sere te variables  could  also  be  introduced 
in  the  problem  statement. 


(2.2) 

(2.3) 
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The  structural  optimization  problem  embodied  in  Eqs.  (2.1- 
2.3)  is  a  nonlinear  mathematical  programming  problem  to  which 
standard  minimization  techniques  can  be  applied.  However  this 
problem  exhibits  some  characteristics  that  make  it  complicated 
when  practical  structural  design  applications  are  considered. 

The  main  difficulty  arises  from  the  fact  that  the  behavior  cons¬ 
traints  (2.2)  are  in  general  implicit  functions  of  the  design 
variables  and  their  precise  numerical  evaluation  for  a  particular 
design  requires  a  complete  finite  element  analysis.  Since  the 
solution  scheme  is  essentially  iterative,  it  involves  a  large 
number  of  structural  reanalyses.  Therefore  the  computational 
cost  often  becomes  prohibitive  when  large  structural  systems 
are  dealt  with. 

The  well  known  necessary  KUHN-TUCKER  conditions  characte¬ 
rize  any  local  minimum  of  the  nonlinear  programming  problem 
(2.  1-2.3)  [  3  ]  .  They  are  written  as  follows  for  each  design 
variable  : 


u  "  _ 1  —  V  +  £  -  Q 

3a.  i  =  1  i  9a,  i  i 

i  j  J  i 


i  =  1  ,  n 


(2.1*) 


with  the  complementary  conditions  : 


r .h .  =0 

J  J 

r  . 

J 

> 

0 

II 

3 

(2.5) 

v.(a.  -  a . ) =  0 

ii  —l 

v  . 

l 

> 

0 

i  =  1  ,  n 

(2.6) 

t.(a.  -  a.)  =  0 

ii  i 

t . 
l 

> 

0 

i  =  1  ,  n 

(2.7) 

The  quantities  (r.,  j  =  1 ,m) ,  associated  with  the  behavioral 
constraints  (2.2),  and  ( v  ^  ,  t  ^  ,  i  =  l,n),  associated  with  the 
side  constraints  (2.3),  are  called  dual  variables  as  opposed 
to  the  primal  variables,  which  are  the  a^.  They  have  the 
meaning  of  lagrangian  multipliers  conjugated  to  the  constraints. 
Depending  upon  whether  a  given  constraint  becomes  an  equality 
or  not  at  the  optimum  (i.e.  is  active  or  inactive),  the  corres¬ 
ponding  dual  variable  is  positive  or  equal  to  zero  (see  Eqs. 
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2.S-2.7).  The  KlIHN-TUCKER  conditions  are  in  general  only 

necessary.  In  the  special  case  of  a  convex  problem,  they 

become  also  sufficient  and  characterize  a  global  optimum. 

They  ran  then  be  used  to  relate  the  primal  variables  -  i.e. 

design  variables  a.  -  to  the  dual  variables  -  i.e.  lagrangian 

multipliers  r.,  v.,  t.  - 
J  i  i 


In  the  last  decade  essentially  two  main  approaches  have 
been  used  to  solve  the  problem  (2. 1-2.3).  One  is  based  on 
the  many  rigorous  numerical  methods  of  nonlinear  mathematical 
programming.  The  other  uses  the  more  intuitive  concepts  of 
optimality  criteria.  These  approaches  have  often  been  opposed 
in  the  past  and  two  corresponding  schools  developed.  The 
advantages  claimed  for  the  mathematical  programming  methods 
are  their  sound  foundations,  their  convergence  properties 
which  can  most  often  be  guarantied  and  their  generality  which 
permits  consideration  of  any  type  of  constraints.  Their  es¬ 
sential  disadvantage  lies  in  the  computing  time  which  increases 
rapidly  with  the  size  of  the  problem,  leading  to  unacceptable 
cost  even  for  relatively  simple  problems. 


The  optimality  criteria  are  based  on  explicit  approximations 
of  the  behavior  constraints  which  are  exact  in  special  cases, 
most  of  the  time  in  statically  determinate  cases.  These  appro¬ 
ximations  are  supposed,  on  intuitive  basis,  to  hold  in  the 
general  case.  Then  the  statement  of  the  problem  (2. 1-2.3)  often 
disappears  since  redesign  formulae  can  be  derived  by  various 
means,  like,  for  instance,  the  KUHN-TUCKER  conditions  (2.4-2. 7) 
in  which  the  explicit  approximations  are  introduced.  The  ad¬ 
vantages  of  using  the  optimality  criteria  are  the  computing  cost 
which  is  low  and  the  intuitive  basis,  which  is  often  appealing 
to  the  engineer.  In  fact  the  number  of  reanalysis  cycles  does 
not  increase  with  the  complexity  of  the  structure,  nor  with  the 
number  of  design  variables.  As  a  consequence  of  this  property, 
untii  recently,  large  scale  examples  of  structural  optimization 
have  only  been  solved  by  optimality  criteria.  The  essential 
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d  i  sadvantages  of  the  optimality  criteria  approach  are  its  lack 
of  generality  and  of  sound  mathematical  foundations,  which  ex¬ 
plains  the  often  impredictable  convergence  properties  and  even 
the  convergence  to  non  optimal  solutions. 

It  has  been  shown  in  a  recent  work  [1,  2,  41  that  these 
two  approaches  are  in  fact  related  to  each  other  and  that  the 
optimality  criteria  can  be  viewed  as  special  mathematical  pro¬ 
gramming  methods  using  linearization  for  the  constraint  surfaces. 
The  main  results  of  the  work  reported  in  Ref.  [  I  1  will  be  sum¬ 
marized  in  the  next  three  sections  of  this  chapter. 

2.2.  The  explicit  approximate  problem 

The  non  linear  programming  problem  to  be  solved  involves 
minimization  of  an  explicit  objective  function  subject  to 
explicit  side  constraints  and  implicit  behavior  constraints. 

As  previously  mentioned  the  implicit  character  of  the  func¬ 
tions  hj  appearing  in  Eq  .  (2.2)  constitutes  the  main  difficulty 

in  achieving  efficiency  in  structural  optimization  methods. 

That  is  why  recoursing  to  explicit  approximation  of  the  behavior 
constraints  appears  to  be  the  only  possible  way  of  successfully 
attacking  the  problem  under  consideration. 


Most  of  the  optimality  criteria  approaches  ileal  with 
problems  involving  constraints  on  static  stresses  and  dis¬ 
placements,  in  which  case  the  behavior  constraints  of  Kq  . 
('?..?)  can  he  written 


h.(a)  =  u.  -  u.(a)  >  0 
1  Id 


(2.  ft) 


where  denotes  an  upper  bound  to  a  response  quantity  u^(a) 
(stress,  nodal  displacement,  relative  displacement,  etc...). 
I' sing  virtual  load  considerations,  explicit  approximations 
of  the  behavior  constraints  (2.8)  can  be  generated  : 


j  (  a  ) 


iE=, 


J_! 

a  . 


>  0 


(2.9) 


where  the  coefficients  c„  are  related  to  virtual  energy  den¬ 
sities  in  the  structural  members. 


The  way  theses  coefficients  are  computed  is  explained 

in  detail  in  Ref.  f  1  )  and  will  be  briefly  reviewed  in  chapter 

3.  For  the  moment  it  is  sufficient  to  recall  that  the  c..'s 

ij 

are  constant  coefficients  in  the  case  of  a  statically  determi¬ 
nate  structure,  so  that  Eq.  (2.9)  represents  then  the  exact 
explicit  form  of  the  behavior  constraints.  In  the  case  of  a 
statically  indeterminate  structure,  the  c^'s  depend  implicitly 
on  the  design  variables  because  structural  redundancy  leads  to 
redistribution  of  the  internal  forces  when  the  member  sizes  are 
modified.  Therefore  the  explicit  constraints  embodied  in  F.q  . 
(2.9)  constitute  in  general  approximate  forms  of  the  original 
constraints  (2.8).  As  shown  in  Ref.  |  1  ]  the  basic  idea  in  the 
optimality  criteria  approach  can  be  viewed  as  transforming  the 
initial  implicit  problem  into  a  sequence  of  explicit  problems. 
Each  explicit  problem  results  from  replacing  the  behavior  cons¬ 
traints  (2.8)  bv  their  approximate  expressions(2.9) . 

On  the  other  hand  the  mathematical  programming  approach, 
after  a  period  of  unefficiency,  has  finally  evolved  into  a 
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powerful  and  now  well  established  design  procedure  that  is  also 
based  upon  explicit  approximation  of  the  behavior  constraints. 

Such  an  approach  was  developed  independently  by  SCHMIT  [  5,6,7,81 
as  the  "approximation  concepts  approach"  and  by  FLEllRY  [  1, 9,  10,  111 
as  the  "mixed  method".  Both  methods  proceed  by  constructing  a 
high  quality  explicit  approximation  of  the  initial  problem  and 
solving  it  partially,  using  a  primal  mathematical  programming 
algorithm,  before  reanalyzing  the  structure  and  updating  the  ap¬ 
proximate  problem  statement.  This  process  facilitates  generation 
of  a  sequence  of  steadily  improved  feasible  designs,  which  is  an 
attractive  feature  for  practical  design  purposes. 


The  key  idea  in  the  mixed  method  [  1,9-11]  ,  as  well  as  in  the 
approximation  concepts  approach  [  5-8] ,  is  to  linearize  the  beha¬ 
vior  constraints  with  respect  to  the  reciprocal  design  variables 


(2.  10) 


Justification  for  this  change  of  variables  lies  in  the  fact  that 
the  stresses  and  the  displacements  are  linear  functions  of  the 
x^'s  for  a  statically  determinate  structure.  Therefore  in  a  mo¬ 
derately  hyperstatic  case,  it  can  be  expected  that  the  constraint 
surfaces  are  very  shallow  and  close  to  planes  in  the  reciprocal 
design  variable  space.  It  is  then  possible,  in  a  method  of  pro¬ 
jected  gradient,  to  progress  with  much  larger  steps  along  tangent 
planes  without  seriously  violating  the  constraints.  Moreover  the 
shallowness  of  the  constraint  surfaces  implies  that  their  linea¬ 
rized  forms  are  usually  very  good  approximations  and  therefore 
the  structure  does  not  need  to  be  reanalyzed  after  each  iteration 
in  a  mathematical  programming  algorithm.  The  linearized  behavior 
constraints  are  obtained  using  first  order  Taylor  series  expansion 
in  terms  of  the  reciprocal  variables  x.  : 


ft.(x) 


n  d  u  . 

-  (  u?  +  .Z  (x .  ~  x?) 

1=1  9  x  .  1  1 

J  1 


0  (2.11) 


where  the  upperscript  °  denotes  quantities  evaluated  at  the 
actual  design  point  x°,  where  the  structural  analysis  is  per¬ 
formed.  Note  that  the  finite  element  analysis  capability 
must  include  auxiliary  sensitivity  analyses  for  evaluating 
first  partial  derivatives  of  the  response  quantities.  Most 
often  the  well  known  pseudo-loads  technique  is  employed  (see 
section  3  . 1. 2 )  . 


Now  it  can  be  shown  f  see  Refs.  I  and  4  ]  that  the  explicit 
approximations  of  the  behavior  constraints  used  in  both  the  op¬ 
timality  criteria  and  mathematical  programming  approaches  (Eqs. 
2.9  and  2.11,  respectively)  are  identical.  Indeed  the  virtual 
energy  densities  c..  employed  in  the  optimality  criteria  appro¬ 
aches  are  nothing  else  that  the  gradients  of  the  response  quan¬ 
tities  with  respect  to  the  reciprocal  variables  : 


c  .  . 
i.l 


au . 

_ I 

3x . 

l 


(2.12) 


Furthermore  the  definition  of  the  c^j’s  following  from  the  vir¬ 
tual  load  technique  clearly  indicates  that 


o  _  ?  o  o 

U  .  -  .  I  .  C  .  .  X  . 

J  1=1  1J  1 


(2.13) 


Therefore  F,q  .  (2.11)  can  be  rewritten 


%  —  ^  O 

h  .  (x)  i  u  .  -  Y  c  .  .  x .  >  0 
J  1  i  =  1  i  J  i 


(2.14) 


which  is  equivalent  to  F. q  .  (2.9). 

In  conclusion  a  unified  structural  optimization  approach  has 
emerged,  which  consists  in  replacing  the  initial  problem  (2. 1-2.3) 
with  a  sequence  of  explicit  approximate  -  or  linearized  -  problems 
Each  linearized  problem  exhibits  the  following,  form  when  written 
in  terms  of  the  reciprocal  variables  : 


minimize 


W 


(2.15) 
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n 


i£i 


p.  i . 


i _ i 

X  . 


1 


subject  to 


n 

u  .  -  .  Z 

J  i* 

.  c  .  .  x  .  >0 

1  1  J  * 

j  =  1  ,m 

(2. 16) 

x  .  <  x  . 

—  l  l 

<  X. 

1 

i  =  1  ,  n 

(2.17) 

where  x.  =  1/a.  and  x.  =  1/a.  are  the  new  side  constraints. 

- 1  i  l  —l 

Tt  is  important  to  mention  that  this  basic  approach  of  trans¬ 
forming  the  initial  problem  into  a  sequence  of  explicit  problems 
is  now  widely  recognized  and  it  is  routinely  employed  for  large 
scale  industrial  applications  (see  PETIAU  [12  ]). 

Various  solution  schemesare  available  to  treat  the  explicit 
problem  embodied  in  Eqs.  (2.15-2.17).  Depending  upon  their  pri¬ 
mal  or  dual  character  the  convergence  properties  of  the  whole  op¬ 
timization  process  will  be  different. 


2.3.  The  mixed  method 


The  linearized  problem  (2.15-2.17)  is  still  a  nonlinear 
programming  problem  (because  the  objective  function  keeps  its 
non  linear  character),  however  it  is  now  explicit  and  easily 
treated  by  standard  minimization  techniques.  In  order  to 
maintain  a  primal  philosophy  (sequence  of  steadily  improved 
feasible  designs),  the  approximation  concepts  approach,  as 
initially  proposed  by  SCHMIT  and  MIURA  [  6,7  ]  ,  employed  either 
a  feasible  direction  method  or  an  interior  penalty  function 
method  to  solve  partially  the  linearized  problem.  In  this  way 
it  is  possible  to  preserve,  at  each  stage  of  the  process,  the 
feasibility  of  the  design  point  with  respect  to  the  primary 
problem  (2.  1-2.3).  Later  the  interior  penalty  function  method 
(called  NEWSUMT)  only  was  retained  [8,13  ].  So  the  NF.WSUMT  op¬ 
timizer  of  the  ACCESS-3  program  f  13,14  ]  tends  to  generate  a 
sequence  of  design  points  that  "funnels  down  the  middle"  of  the 
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feasible  region.  This  represents  an  attractive  feature  from 
an  engineering  point  of  view  since  it  provides  the  design  with 
acceptable  non  critical  designs  at  each  stage  of  the  process. 
Furthemore  the  convergence  properties  of  the  method  can  be  con¬ 
trolled  by  modifying  the  number  of  response  surfaces  (typically 
1  or  2)  and  the  response  factor  decrease  ratio  (typically  0.5 
through  0.1). 

On  the  other  hand,  starting  from  an  optimality  criteria 
approach,  a  method  similar  to  the  approximation  concepts  ap¬ 
proach  was  independently  developed  in  Refs.  [9-11  ].  This 
mixed  method  uses  virtual  load  considerations  to  generate  a 
first  order  approximate  problem  which  is  identical  to  the  li¬ 
nearized  problem  posed  by  Eqs.  (2.15-2.17).  This  problem  is 
also  solved  partially  using  a  primal  solution  scheme,  with  the 
aim  of  preserving  the  design  feasibility  as  in  the  approximation 
concepts  approach.  Each  iteration  in  the  mixed  method  is  made 
up  of  a  restitution  phase  and  a  minimization  phase. 

The  restitution  phase  consists  in  bringing  the  design  point 
back  to  the  boundary  of  the  feasible  region.  It  is  based  upon 
the  concept  of  scaling,  which  is  classical  in  the  optimality 
criteria  approaches.  Scaling  does  not  introduce  any  redistri¬ 
bution  of  the  internal  forces  in  the  structure  and  therefore 
it  permits  generation  of  a  feasible  boundary  point  without  any 
additional  structural  analysis.  In  the  design  space  scaling 
corresponds  to  a  move  along  a  straight  line  joining  the  origin 
to  the  current  design  point,  where  the  structural  analysis  is 
made  (sec  Fig.  I).  Along  the  scaling  line  the  coefficients  c  ^  j 
(see  F.q  .  2.9)  remain  constant  and  the  explicit  approximations 
lead  to  the  exact  values  of  the  constraints  and  their  gradients 
(see  section  3.1.1).  As  a  result  the  linearized  problem  statement 
is  not  affected  by  scaling.  Geometrically  each  real  restraint 
surface  is  approximated,  in  the  space  of  the  reciprocal  variables, 
by  its  tangent  plane  at  its  point  of  intersection  with  the  sca¬ 
ling,  line  (see  Fig.  1 ) . 
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Starting  from  a  feasible  boundary  point,  the  minimization 
phase  in  the  mixed  method  consists  in  reducing  the  objective 
function  in  the  subspace  defined  by  the  tangent  planes  to  the 
constraint  surfaces.  This  is  achieved  by  solving  the  linea¬ 
rized  problem  (2.15-2.17)  using  a  projection  algorithm,  but 
stopping  the  minimization  after  a  limited  number  of  steps, 
denoted  k,  before  reaching  the  minimum  of  the  linearized  pro¬ 
blem.  After  ic  one -d imens ional  minimizations  have  been  ac¬ 
complished,  the  structure  is  reanalyzed,  scaling  is  performed, 
the  linearized  problem  is  reformed  and  again  solved  partially. 

That  optimization  process  is  illustrated  on  Fig.  2  in  a 
hypothetical  2-dimensional  space.  When  k  is  limited  to  one 
stop,  the  method  is  a  strict  application  of  a  gradient  pro¬ 
jection  type  of  algorithm  to  the  original  problem  (2.  1-2.3). 
The  method  then  exhibits  all  the  properties  of  the  primal 
mathematical  programming  approaches,  that  is  high  cost  but 
guarantied  convergence.  When  k  is  not  limited,  the  linea¬ 
rized  problem  is  solved  completely  before  reanalyzing  the 
structure  and  updating  the  linearized  problem  statement.  The 
method  then  becomes  an  optimality  criteria  approach  (see  next 
section).  The  corresponding  properties  of  fast  but  uncertain 
convergence  are  to  be  expected.  Finally,  when  k  is  limited 
to  a  given  finite  number,  the  method  is  mixed.  This  leads  to 
considering  the  number  of  steps  k  as  a  convergence  control  pa¬ 
rameter  that  should  be  assigned  high  values  for  economy  and  re¬ 
duced  when  divergence  occurs. 

As  implemented  in  the  SAMCEF  program  (  15,16  ] ,  the  minimi¬ 
zation  algorithms  used  in  the  mixed  method  are  based  either  on 
the  conjugate  gradient  method  with  an  orthogonal  projection  o- 
perator  or  on  the  generalized  Newton's  method  with  a  weighed 
projection  operator.  This  latter  method  has  also  been  intro¬ 
duced  in  the  ACCF.SS-3  program  [  13,14  I  . 
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From  a  mathematical  programming  point  of  view,  the  appro¬ 
ximation  concepts  approach  and  the  mixed  method  can  bo  clas¬ 
sified  as  mixed  primal-linearization  methods.  The  initial 
problem  is  transformed  into  a  sequence  of  linearized  problems, 
which  is  classical  in  the  linearization  methods  of  mathematical 
programming.  However  each  subproblem  is  solved  only  partially 
us-ing  a  primal  solution  scheme  that  insures  feasibility  of  the 
intermediate  designs  at  each  stage. 


2.4.  The  generalized  optimality  criteri o n 

Tn  the  previous  section  a  primal  philosophy  has  been  adopted 
that  leads  to  partial  solution  of  the  linearized  problem  (2.15- 
2.17)  using  an  interior  penalty  function  formulation  with  a 
small  number  of  response  surfaces  or  a  projection  algorithm  with 
a  small  number  of  one-dimensional  minimization.  This  primal  so¬ 
lution  scheme  produces  a  sequence  of  feasible  designs  with  de¬ 
cree  s  i  ng  va  1  ue  s  of  the  objective  function.  An  alternative  view¬ 
point  is  to  recognize  that  the  approximation  made  by  linearizing 
the  constraints  with  respect  to  the  reciprocal  design  variables 
is  of  such  high  quality  that  the  current  explicit  problem  can 
be  solved  exactly,  rather  than  partially,  after  each  structural 
ronnnlysis.  This  idea  leads  to  abandoning  the  primal  philosophy 
in  favor  of  a  pure  linearization  approach  in  mathematical  pro¬ 
graming,  but  it  can  also  be  viewed  as  adopting  an  optimality  cri¬ 
teria  strategy. 


2.4. I .  Dual  solution  scheme 

Since  only  the  final  exact  solution  of  each  linearized  pro¬ 
blem  needs  to  be  known  at  each  redesign  stage,  any  minimization 
algorithm  can  be  chosen  to  solve  it.  In  order  to  improve  the 
nmput.iti  onal  efficiency,  it  is  advisable  to  select  a  specialized 
nonlinear  programming  algorithm,  well  suited  to  the  particular 
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m;i  l  h  >  •  in  .'i  t  i  c  n  I  s  l  rue  hire  of  the  explicit  problem  (  2  .  I  S  - 2 .  I  7  )  . 
The  objective  function  is  strictly  convex  and  the  constraints 
are  linear,  so  that  the  problem  (15-17)  is  a  convex  program¬ 
ming  problem.  Moreover  all  the  functions  defining  this  pro¬ 
blem  are  explicit  and  separable.  In  such  a  case  at  least 
two  kinds  of  method  seem  to  be  efficient  :  the  second  order 
primal  projection  methods  and  the  dual  methods!  see  Kef.  17  ]  . 


The  dual  method  formulation  is  especially  attractive,  be¬ 
cause  the  dual  problem  presents  a  much  simpler  form  than  the 
primal  problem.  Indeed,  for  a  convex  problem,  the  lagrangian 
multipliers  associated  with  the  constraints  have  the  meaning 
of  dual  variables  in  terms  of  which  an  auxiliary  and  equivalent 
problem  can  be  stated.  Under  some  unrestrictive  conditions, 
this  dual  problem  can  be  reduced  to  the  maximization  of  the  la¬ 
grangian  functional  with  simple  nonnegativity  requirements  on 
the  dual  variables.  Since,  in  addition,  the  explicit  problem 
(2.15-2.17)  is  of  separable  form,  the  dual  formulation  leads  to 
a  very  efficient  solution  scheme  since  each  primal  variable  can 
be  analytically  expressed  in  terms  of  the  dual  variables.  Futher- 
more  the  dimensionality  of  the  dual  problem  is  equal  to  the  num¬ 
ber  of  linearized  behavior  constraints  (2.16),  which  is  most 
often  small  when  compared  to  the  number  of  design  variables. 
Therefore  the  dual  problem  exhibits  a  simpler  form  and  a  lower 
dimensionality  than  the  primal  problem. 

From  a  mathematical  programming  point  of  view,  the  dual 
methods  are  thus  likely  to  provide  the  most  efficient  solution 
scheme  to  the  linearized  problem  (2.15-2.17),  provided  this 
problem  may  be  solved  exactly  at  each  stage,  rather  than  par¬ 
tially,  using  a  primal  method  [see  Ref.  18  ].  This  demands  that 
the  original  behavior  constraints  exhibit  a  moderate  nonlinea¬ 
rity  with  respect  to  the  reciprocal  variables,  which  is  actually 
true  for  problems  involving  stress,  displacement,  frequency  and 
buckling  constraints.  Another  significant  advantage  of  the  dual 
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methods  is  that  they  can  take  into  c ons i der a t i on  discrete  design 
variables  without  weakening  the  efficiency  of  the  optimization 
process  [13,19  ] .  Discrete  variables  are  useful  for  representing 
comercially  available  gage  sizes  of  sheet  metal,  the  number  of 
plies  in  a  laminated  composite  skin,  etc...  Finally  duality 
techniques  can  also  be  employed  to  monitor  the  progress  of  the 
optimization  process,  as  proposed  by  MORRIS  et  al.  [20  ]  . 

Two  maximization  algorithms  are  available  in  the  SAMCEF  pro¬ 
gram  to  solve  the  dual  problem  :  a  first  order  conjugate  gradient 
type  of  algorithm  and  a  second  order  Newton  type  of  algorithm 
[1,16  ]  .  Recently  these  dual  algorithms  were  implemented  in  the 
ACCESS-3  program,  which  is  based  on  the  approximation  concepts 
approach  [  13,14  ] .  The  second  order  algorithm  was  improved  by 
employing  a  simple  but  efficient  one-dimensional  maximization 
scheme.  The  first  order  one  was  deeply  modified  and  it  is  now 
capable  of  treating  problems  involving  discrete  design  variables. 


2.4.2.  Explicit  redesign  relations 

Turning  now  to  an  optimality  criteria  strategy,  the  simple 
form  of  the  linearized  problem  (see  Eqs.  2.15-2.17)  suggeststo 
make  use  of  its  explicit  character  in  order  to  derive  redesign 
formulas  for  the  design  variables.  In  fact,  as  shown  in  Ref. 

I  4  1  ,  the  whole  process  of  combining  the  linearization  of  the 
behavior  constraints  with  respect  to  the  reciprocal  design 
variables  and  a  dual  solution  scheme  can  be  viewed  as  a  genera¬ 
lization  of  the  optimality  criteria  approaches.  This  process 
consists  in  solving  exactly,  after  each  structural  reanalysis, 
the  linearized  problem  (2.15-2.17),  which  can  be  recast  as  fol¬ 
lows  in  terms  of  the  direct  design  variables  : 


n 


minimi z  p 

W  =  .  T. ,  p  .  1  .  a  . 

l  =  1  l  l  l 

(2.18) 

subject  to 

n  c  i  ; 

u.  -  .  Z  .  — >  0 

j  * 

1  ,m 

(2.19) 

J  i=l  a. 

a  .  >  a  .  >  a  . 
l  l—i 

i  * 

1  ,  n 

(2.20) 
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Instead  of  employing  primal  or  dual  mathematical  programming 
methods,  an  alternative  approach,  which  is  typical  of  the  opti¬ 
mality  criteria  philosophy,  is  to  use  the  explicit  character  of 
the  approximate  problem  (2.18-2.20)  in  order  to  express  analy¬ 
tically  the  optimal  design  variables.  This  can  be  achieved  using 
the  KUHN-TUCKER  conditions  (2.4-2. 7)  which,  in  view  of  the  con¬ 
vexity  of  the  problem,  are  sufficient  for  global  optimality. 

These  conditions  lead  to  a  generalized  optimality  criterion 
yielding  explicitly  the  design  variables  : 


active  design  variables  : 

1  1/2  2  , 

a .  =( -  E  c  .  .  r.)  if  0.Jl.a.  <"  £  c  .  .  r.  <  n.fc.af 
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p assive  design  variables  : 


a . 

if 

£ 

c  .  . 
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a  . 
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E 
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1  J 

J 
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(2.21) 


(2.22) 

(2.23) 


In  these  expressions,  the  lagrangian  multipliers  r.  are  as¬ 
sociated  with  the  linearized  behavior  constraints  (2.19).  They 
must  satisfy  the  complementarity  conditions  (2.5),  namely  : 


active  constraint  : 


r .  >  0  if 
J 


inactive  constraint  : 


r .  =  0  if 
J 


c  .  . 


(2.24) 


(2.25) 


The  equations  (2.21-2.23)  relating  the  design  variables  a  - 
to  the  lagrangian  multipliers  r^  provide  a  basis  for  separating 
the  design  variables  in  two  groups.  The  passive  variables  are 
those  that  are  fixed  to  a  lower  or  an  upper  limit  (see  F.qs.  2.22 
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and  2.23,  respectively)  while  the  active  variables  are  explici¬ 
tly  given  in  terms  of  the  lagrangian  multipliers  using  Eq .  (2.21). 

This  subdivision  of  -he  design  variables  into  active  and  passive 
groups  is  classical  in  the  optimality  criteria  approaches  f  21-29] . 
Equations  (2.22,  2.23)  clearly  indicate  that  the  dual  space  - 
i.e.  the  space  of  the  lagrangian  multipliers  r^  -  is  partitioned 
into  several  regions  corresponding  to  different  divisions  in 
active  and  passive  design  variables.  These  regions  are  separated 
from  each  other  by  planes  across  which  the  second  derivatives 
of  the  lagrangian  function  are  discontinuous  (see  Refs.  [  1,30  ] 
for  more  details)  . 

When  the  lagrangian  multipliers  satisfying  F.qs.  (  2.24,2.25) 
are  known,  the  optimal  design  variables  can  be  easily  computed 
using  the  explicit  optimality  criterion  (2.21-2.23).  Therefore 
the  problem  has  been  replaced  with  a  new  one,  which  is  defined 
in  terms  of  the  lagrangian  multipliers  only.  To  solve  this  new 
problem,  the  conventional  optimality  criteria  techniques  usually 
make  the  assumption  that  the  set  of  active  constraints  is  known 
in  advance,  avoiding  thus  the  inequality  constraints  on  the  la¬ 
grangian  multipliers  appearing  in  Eqs.  (2.24,2.25).  An  update 
procedure  for  the  retained  lagrangian  multipliers  is  then  em¬ 
ployed,  so  that  the  optimal  design  variables  can  be  sought  i- 
teratively  by  coupling  the  update  procedure  and  the  explicit 
optimality  criterion  (2.21-2.23).  The  essential  difficulties 
involved  in  applying  these  optimality  criteria  methods  are  those 
associated  with  identifying  the  correct  set  of  active  constraints 
and  the  proper  corresponding  set  of  passive  members  [  25,28  ] . 

These  difficulties  were  recognized  and  addressed  with  varying 
degrees  of  success  in  studies  such  as  those  reported  in  Refs, 
f  26,27  ] .  However  it  was  only  with  the  dual  formulation  set 
forth  in  Refs.  [  1,30  ]  that  these  obstacles  were  conclusively 
overcome . 

The  dual  method  approach  consists  in  maximizing  the  lagrangian 
function  subject  to  nonnegativity  constraints  on  the  lagrangian 
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multipliers.  This  approach  can  therefore  be  viewed  as  using 
an  update  procedure  for  the  lagrangian  multipliers,  exactly 
like  the  conventional  optimality  criteria  techniques.  After 
the  update  procedure  is  completed,  the  design  variables  can 
be  evaluated  using  the  optimality  criteria  equations  (2.21- 
2.23).  Since  the  dual  maximization  problem  exhibits  a  re¬ 
markably  simple  form,  its  exact  solution  can  be  generated  at 
a  low  computational  cost,  which  is  comparable  to  that  required 
by  the  recursive  techniques  of  conventional  optimality  criteria 
approaches.  The  dual  algorithms  can  handle  a  large  number  of 
inequality  constraints  and  they  intrinsically  contain  a  rational 
scheme  for  identifying  the  active  constraints  through  the  non¬ 
negativity  constraints  on  the  lagrangian  multipliers  (or  dual 
variables).  They  also  automatically  sort  out  the  active  and 
passive  design  variable  groups  using  the  explicit  relation¬ 
ships  between  primal  and  dual  variables. 


When  the  design  variables,  instead  of  varying  continuously, 
can  only  take  on  available  discrete  values,  it  can  be  shown 
that  the  optimality  criteria  equations  (2.21-2.23)  must  read 
as  follows  for  each  discrete  variable  (see  Refs.  [  13]  and  [  19  ]  )  : 


k  k-1 


k  k+  1 


if  p  •  2.  .  a  V  a  V  <  £  c  .  .  r  .  <p  .  I .  aTa? 

p  i  i  i  l  j  ijjVM1iii 


(2.26) 


whore  it  is  understood  that  {a.,  k  =  1,2,...)  denotes  the  set 

of  available  discrete  values  for  the  i  ^  design  variable.  These 
expressions  show  that  the  dual  space  is  subdivided  into  several 
regions,  each  of  which  corresponds  to  a  distinct  combination 
of  available  discrete  values  of  the  design  variables.  These 
regions  are  separated  from  each  other  by  planes  across  which 
the  first  derivatives  of  the  dual  function  are  discontinuous 
(see  Ref.  1  13,19  ]  for  more  details). 


Tn  summary  a  generalized  optimality  criteria  approach  can 
be  defined  in  the  mathematical  programming  terminology  as  a 
special  form  of  the  linearization  methods.  It  amounts  to  re- 
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placing  the  original  problem  with  a  sequence  of  linearized 
problems  and  applying  a  dual  solution  scheme  to  each  sub¬ 
problem.  It  should  be  noted  however  that  only  the  behavior 
constraints  are  linearized  with  respect  to  the  reciprocal 
design  variables,  while  in  classical  linearization  methods, 
the  objective  function  is  linearized  too.  This  is  not  ne¬ 
cessary  in  the  present  approach  because  the  structural  weight 
is  an  exact  explicit  function  of  the  reciprocal  variables. 

It  should  be  also  emphasized  that  unlike  the  primal  algorithms 
discussed  in  section  2.3,  the  dual  algorithms  cannot  be  used 
to  solve  only  partially  the  approximate  problem  embodied  in 
Eqs.  (2.18-2.19),  because  intermediate  points  in  the  dual  space 
usually  correspond  to  higly  infeasible  points  in  the  primal 
space.  Therefore  the  possibility  of  controlling  the  convergence 
of  the  optimization  process  disappears  when  a  dual  optimizer  is 
selected  in  the  SAMCEF  [ 15,16  ]  or  ACCESS-3  [ 13,14]  programs. 
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3.  STATIC  CONSTRAINTS  ON  STRESSES  AND  DISPLACEMENTS 


This  chapter  is  concerned  with  structural  optimization 
problems  involving  behavior  constraints  on  static  stresses 
and  displacements.  Most  of  the  material  presented  recalls 
conclusions  drawn  in  Ref.  [  1  ]  ,  however,  some  new  develop¬ 
ments  will  also  be  reported  with  regards  to  applied  loads 
that  depend  on  the  design  variables. 

3.  I  .  Evaluation  of  the  constraint  gradients 


As  emphasized  in  chapter  2,  the  mixed  method  and  the 
generalized  optimality  criterion  set  forth  in  Ref.  [  I  ] 
are  both  based  upon  linearization  of  the  stress  and  dis¬ 
placement  constraints.  Such  a  linearization  process  re¬ 
quires  gradient  evaluations.  Two  different  procedures  are 
available  to  this  end,  namely  the  virtual  load  and  the  pseudo¬ 
loads  techniques. 

3.1.1.  The  virtual  load  technique 

Most  of  the  conventional  optimality  criteria  approaches, 
as  well  as  the  generalized  optimality  criterion  reported  in  ref.^l  j  , 
use  the  virtual  load  technique  to  generate  explicit  approxi¬ 
mations  of  the  stress  and  displacement  constraints.  In  Ref. 
Mia  flexibility  is  defined  as  a  linear  combination  of  dis¬ 
placements  : 


u  . 
J 


(3.1) 


where  b  .  is 
J 

load  vector 
f  1  e  x  i  h  i  l  i  t  y 
sum  of  the 


a  vector  of  constants.  Considering  a  virtual 

given  numerically  by  b.,  it  follows  that  the  j*"'1 
t  h  ^ 

under  the  l  loading  case  can  be  expressed  as  the 
contributions  of  each  element  : 


u 


j* 


(3.2) 


...... 


_ _ _ 
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with 


c . . .  =  (q.  K.  q  )  a. 
1  j  J,  l 


(3.3) 


In  these  expressions  q  ^  and  qf  are  respectively  the  virtual 
and  real  displacement  vectorsand  denotes  the  stiffness 
matrix  of  the  i ^  finite  element  (see  Eq.  1.9).  The  flexi¬ 
bility  coefficients  are  constant  for  a  statically  de¬ 

terminate  structure,  in  which  case  Eq.  (3.2)  furnishes  the 
exact  explicit  expression  of  the  flexibilities.  Therefore 
explicit  optimality  criteria  can  be  derived,  which  define 
analytically  the  optimal  values  of  the  design  variables  when 
displacement  constraints  are  imposed. 


In  a  finite  element  model  the  stresses  are  also  linear 
comb inat ions  of  the  displacements,  just  as  the  flexibilities 
defined  in  Eq .  (3.1).  Indeed  a  given  stress  component  in  the 

k*"*1  element  is  written  as 


a 


k 


(3.4) 


where  t^  is  a  vector  of  constants  corresponding  to  a  line  of  the 
stress  matrix  (see  Eq .  1.8).  Following  out  the  above  pro¬ 

cedure  the  vector  t,  is  introduced  as  a  virtual  loading  case, 

k 

which  permits  generation  of  an  explicit  expression  for  the 
corresponding  stress  component  under  the  load  condition  : 


w  i  th 


n 


i^l 


(3.5) 


dik*  =  (qk  Ki  V  3i 


(3.6) 


[n  these  equations,  q  ^  denotes  the  virtual  displacement  vector 
due  to  the  virtual  load  case  given  by  t^.In  a  plate  element  the 
constraint  is  usually  placed  on  an  equivalent  stress  whose  square 
is  a  quadratic  form  of  the  displacements.  It  has  been  shown  in 
Ref.  |  1  ]  that  an  explicit  expression  of  the  form  given  by  Eq. 
(3.5)  can  still  be  obtained  by  using  special  virtual  load  cases. 
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So  the  explicit  expressions  of  the  stress  and  displacement  cons¬ 
traints  exhibit  now  the  same  form  (see  Eqs  .  3.2  and  3.5). 

Writing  then  under  the  common  notations  stated  in  Fq  .  (2.9),  it 
appears  that  an  optimality  criteria  approach  proceeds  by  repla¬ 
cing  the  exact  behavior  constraints 


u.  <  u.  (3.7) 

J  J 


by  the  explicit  constraints 


u  . 
J 


(3.8) 


Introducing  virtual  strain  energies  associated  with  each  finite 
e lemen  t , 


e  .  . 
1 J 


(3.9) 


the  generalized  optimality  criterion  reviewed  in  section  2.4 
takes  the  form 


m 

.1  .  r  .  £  .  .  =  c-  (3.10) 

J  =  !  J  1 J 


where  the  e^j's  have  the  meaning  of  virtual  energy  densities 


e  .  . 
1 J 


e  .  . 
i  J 


f>.  Z  .  a  . 
xxi 


(3.11) 


Tn  expression  (3.9),  q.  represents  the  virtual  displacement 

J  .  .  th 

vector  for  a  virtual  load  case  conjugated  to  the  j  behavior 
constraint.  Note  that  the  real  load  case  index  i  is  omitted 
for  sake  of  clarity.  It  should  be  recalled  that  the  optima¬ 
lity  criterion  posed  by  Eq .  (3.10)  applies  only  to  the  active 

design  variables  (see  Eqs.  2.21-2.23).  In  the  special  case 
where  only  one  behavior  constraint  is  assigned  (i.e.,  m  =  1  in 
F.  q .  3.10),  the  optimality  criterion  states  that  in  the  optimal 
structure  the  virtual  strain  energy  density  is  the  same  in  each 
e lernent . 


^ ,  i  pi 
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The  main  result  reported  in  Ref.  [  1  ]  is  that  the  virtual 
load  technique  used  in  the  optimality  criteria  approaches  gene¬ 
rates  first  order  explicit  approximations  of  the  stress  and 
displacement  constraints.  As  previously  mentioned  the  flexi¬ 
bility  coefficients  c..  are  constant  in  a  statically  determi¬ 
nate  structure  and  they  depend  implicitly  on  the  design  varia¬ 
bles  in  a  statically  indeterminate  one.  However  it  is  essential 
to  note  that  they  are  not  affected  by  a  scaling  of  the  design, 
that  is  by  a  multiplication  of  all  the  variables  by  a  common 
factor.  In  the  design  space  such  a  scaling  moves  the  design 
point  along  a  line  joining  the  origin  to  the  current  point. 
Therefore  the  explicit  expression  stated  in  Eq .  (3.8)  yields 

the  exact  value  of  the  flexibility  u  ^  all  along  the  scaling  line. 
Geometrically  this  means  that  the  approximate  constraint  surface 

<\j  - 

Uj  =  u.  passes  through  the  point  of  intersection  of  the  corres¬ 
ponding  exact  constraint  surface  u.  =  u^  with  the  scaling  line 
(see  Fig  .  3  )  . 


In  addition,  it  can  be  shown  that  the  explicit  constraint 
(3.8)  furnishes  also  the  exact  derivatives  of  the  flexibility 
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3u  . 

3  u  . 
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(3.12) 


for  any  design  point  on  the  scaling  line  f  see  Ref.  30  ] .  There¬ 
fore  the  explicit  forms  (3.8)  represent  first  order  approximations 
of  the  flexibility  (or  stress)  constraints  at  all  points  along 
the  scaling  line.  Geometrically  the  approximation  made  in  an  op¬ 
timality  criteria  approach  can  thus  be  viewed  as  replacing  each 
exact  constraint  surface  by  a  tangent  surface  at  its  point  of 
intersection  with  the  scaling  line  (see  Fig.  3  ).  From  Eq.  (3.12) 

it  is  now  apparent  that  the  gradients  of  the  stress  and  displace¬ 
ment  constraints  with  respect  to  the  reciprocal  variables  are  gi¬ 
ven  by  the  coefficients  c..  (see  Eq .  2.12).  Also  the  explicit 
approximate  constraints  (3.8)  when  written  in  terms  of  the  reci¬ 
procal  variables  (see  Eq .  2.16)  can  be  identified  as  being  iden- 
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tical  to  the  first  order  Taylor  series  expansion  embodied  in 
Eq  .  (2.11).  This  supports  the  conclusion  stated  in  section 

2.2,  according  to  which  a  unified  approach  to  structural  op¬ 
timization  consists  in  solving  the  original  problem  embodied 
in  Eqs.  (2. 1-2.3)  as  a  sequence  of  linearized  problems  of  the 
form  posed  in  Eqs.  (2.15-2.17). 


From  the  foregoing  developments,  it  appears  that  the  vir¬ 
tual  load  technique  can  be  considered  as  a  particular  procedure 
for  computing  the  constraint  graidents.  Introducing  virtual 
load  cases  b^  for  each  flexibility  constraint  the  corresponding 
virtual  displacement  vectors  are  evaluated  by  solving  the  sys¬ 
tems  of  linear  equations 


K  q 


b  . 
J 


j  =  I 


(3.13) 


The  flexibility  gradients  are  then  computed  as  follows  : 

3  u 
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(3.  14) 


By  the  same  operation  the  explicit  forms  of  the  flexibility 
constraints  defined  by  Eq .  (3.2)  are  available.  Of  course  the 

same  arguments  hold  for  the  stress  components  expressed  in 
Eq.  (3.4)  and  theit  gradients  are  calculated  by  introducing 
virtual  loading  cases  t^  (see  Eqs.  3.5  and  3.6),  yielding 
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(3.15) 


This  approach  to  the  evaluation  of  the  constraint  gradients 
requires  as  many  additional  virtual  loading  cases  as  the  number 
of  stress  and  displacement  constraints  retained  in  the  linearized 
problem  statement,  regardless  of  the  number  of  design  variables 
and  of  the  number  of  real  loading  conditions. 
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3 . I .  2  .  The  pseudo-loads  technique 

On  the  other  hand  most  of  the  mathematical  programming 
approaches,  such  as  the  approximation  concepts  method  deve¬ 
loped  by  SCHMIT  [  5-8,18  ] ,  use  the  pseudo-loads  technique 
of  FOX  [  31  ]  for  computing  the  gradients  of  the  nodal  dis¬ 
placements  characterizing  the  finite  element  model.  The 
gradients  of  the  stress  and  displacement  constraints  are 
then  readily  evaluated. 

Differentiating  the  equilibrium  equations  for  the  H*"*1 
load  case, 

Mt  1  gt  (3.16) 


leads  to 


(3.17) 


In  the  bar  and  membrane  elements  considered  in  this  chapter 
the  stiffness  matrix  is  linear  in  the  design  variables  so  that 


3a.  a . 

l  l 


(3.18) 


whore  K.  is  the  stiffness  matrix  of  the  i*"*1  element  (see  Eq.1.9). 
i 

Therefore  the  nodal  displacement  gradients  can  be  computed  by 
introducing  additional  loading  cases,  called  pseudo-load  vectors: 


(3.19) 


and  solving  the  systems  of  linear  equations 


i  =  1  ,  n 

i.  =  i  ,  c 


(3.20) 


Recalling  the  definitions  (3.1)  of  the  flexibilities  and  (3.4) 
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of  the  stresses,  the  gradients  of  the  behavior  constraints  arc- 
evaluated  fr  om 


and 


3a  . 
1 


T  ^  ^  i 
1  i  3a. 


3o 


kJl 


3a . 

l 


t  ^  ^  a 
'k  3a. 


(3.2  1) 


(3.22) 


The  number  of  pseudo-load  vectors  is  directly  related  to 
the  number  of  design  variables  and  to  the  number  of  applied 
loading  conditions,  and  it  is  independent  of  the  number  of 
behavior  constraints. 


3.1.3.  Selection  of  a  procedure 

The  decision  as  to  which  procedure  should  be  selected  to 
coupute  the  constraint  gradients  can  be  based  on  a  comparison 
of  the  total  number  of  additional  virtual  or  pseudo-load  cases 
introduced  in  the  structural  reanalysis  at  each  given  stage  : 

-  if  the  virtual  load  technique  is  used,  the  number  of 

additional  loading  cases  is  equal  to  the  number  of 
stress  and  displacement  constraints  retained  in  the 
linearized  problem  statement  for  the  current  stage 
(see  F,q.  3.13)  ; 

-  if  the  pseudo-loads  technique  is  adopted,  the  number 
of  additional  loading  cases  is  equal  to  the  number  of 
independent  design  variables  times  the  number  of  ap¬ 
plied  loading  conditions  (see  Eq .  3.20). 

Tt  should  be  noted  that  the  conventional  optimality  cri¬ 
teria  methods  generally  take  into  account  a  large  number  of 
design  variables  and  a  small  number  of  critical  flexibility  cons¬ 
traints,  in  which  case  the  virtual  load  technique  is  obviously 
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the  most  efficient  scheme  for  generating  their  explicit  appro¬ 
ximate  form  (and  thus  for  computing  their  gradients).  On  the 
other  hand  most  of  the  mathematical  programming  approaches  are 
restricted  to  a  relatively  small  number  of  independent  design 
variables,  but  need  evaluation  of  the  gradients  for  all  the 
stress  and  displacement  constraints,  and  the  pseudo-loads 
technique  is  then  more  advantageous.  In  fact  these  two  distinct 
procedures  for  computing  the  constraint  gradients  constitute 
one  of  the  main  reasons  why  optimality  criteria  and  mathematical 
programming  approaches  have  been  traditionally  opposed  in  the 
past. 


Despite  the  apparent  opposition  between  the  two  procedures, 
they  can  be  analytically  related  to  each  other  in  a  very  simple 
way.  Substituting  Eq .  (3.20)  into  Eq  .  (3.21)  yields  the  flexi¬ 

bility  gradients  resulting  from  the  pseudo-load  technique  : 


pi£ 


(3.23) 


Using  the  definition  (3.19)  of  the  pseudo-load  vectors  it 
follows  that 


l 


K 


Ki(ljl 


(3.24) 


Finally  in  view  of  Eq.  (3.13),  it  is  apparent  that 


1  T 

inr~  =  "  it  qj  Ki^ 

i  i 


(3.25) 


which  is  equivalent  to  the  definition  (3.14)  of  the  flexibility 
gradients  derived  from  the  virtual  load  procedure.  Note  that 
the  same  developments  can  obviously  be  made  for  the  stresses  by 
replacing  ( b . ,  q  .  )  with  ( t  k  ,  q  k  )  . 


3.2.  Conventional,  generalized  and  hybrid  optimality  criteria 


The  generalized  optimality  criterion  (COC)  discussed  in 
section  2.4  as  well  as  the  mixed  method  and  the  approximation 
concep ts approach  (section  2.3),  use  a  linearized  expression 
for  each  behavior  constraint.  However,  in  many  conventional 
optimality  criteria  (COC),  such  as  those  reported  in  Refs. 

[  2  1,  22  ,  24,28,29  ]  ,  only  the  displacement  constraints  are  re¬ 
placed  with  first  order  explicit  approximation  while  the  stress 
constraints  are  treated  using  the  classical  "Fully  Stressed 
Design"  (FSD)  concept.  In  this  approach  the  stress  constraints 
embodied  in  Eq  .  (1.4)  are  transformed  into  simple  side  cons¬ 

traints  : 


(3.26) 


The  minimum  values  a^  are  given  by  the  well  known  stress  ratio 
formula  : 


'V 

a  . 
—  i 


o 

a  .  max 


(3.27) 


where  a?  denote  the  design  variables  at  the  current  stage  and 
o?p  the  corresponding  stresses. 


As  shown  in  Ref.  [31  ]  the  FSD  procedure  can  be  interpreted 
as  using  zero  order  approximation  of  the  stresses  on  the  scaling 
line,  because  it  relies  on  explicit  expressions  that  preserve 
only  the  value  of  the  stresses  along  that  line,  and  not  their 
derivatives.  Geometrically,  the  COC  approximates  the  stress 
constraint  surfaces  by  tangent  surfaces,  while  the  FSD  procedure 
uses  planes  normal  to  the  axes  of  the  design  space.  In  both 
approaches  each  approximate  restraint  surface  passes  by  the  point 
of  intersection  of  the  scaling  line  with  the  corresponding  real 
restraint  surface.  A  graphical  interpretation  of  the  differences 
between  conventional  (COC)  and  generalized  (COC)  optimality  cri¬ 
teria  is  illustrated  in  Fig.  3  .  The  analysis  point  is  denoted  a°. 
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the  solutions  of  the  approximate  problems,  a  ,  and  the  solution 
of  the  r  r  a  1  problem,  a  .  This  one  lies  at  the  intersection  of 
the  two  restraint  surfaces,  namely,  a  stress  limitation  in 
element  2  and  a  flexibility  constraint.  In  the  COC  each  appro¬ 
ximate  restraint  surface  is  tangent  to  the  real  one  at  its  point 
of  intersection  with  the  scaling  line.  In  the  COC  this  feature 
applies  only  to  the  flexibility  constraint  u  <  u,  while  the 
stress  limitation  0£  ^  is  represented  in  the  design  space  by 
a  plane  a^  =  a ^  normal  to  the  a ^  axis. 

From  a  computational  point  of  view,  recoursing  to  the  FSD 
procedure  offers  two  important  advantages.  First,  when  the 
virtual  load  technique  is  employed  to  compute  the  constraint 
gradients,  the  number  of  addtional  loading  cases  is  substantially 
reduced,  because  the  zero  order  approximations  of  the  stress 
constraints  de  not  require  any  virtual  load  cases  (see  Eqs.  3.26, 
3.27).  Secondly  the  number  of  behavior  constraints  retained  in 
each  approximate  problem  statement  (see  Eq .  2.19)  is  also  signi¬ 
ficantly  reduced,  because  all  the  stress  constraints  are  now 
replaced  with  side  constraints.  This  feature  is  especially  be¬ 
neficial  when  a  dual  solution  scheme  is  used  to  solve  the  explicit 
problem,  since  the  dimensionality  of  the  dual  problem  is  equal 
to  the  number  m  of  linearized  constraints  embodied  in  Fq .  (2.19). 

On  the  other  hand,  it  is  well  known  that  the  FSD  procedure  does 
not  always  converge  to  the  true  optimum  and  sometimes  leads  to 
instability  or  divergence  of  the  optimization  process. 

So  we  have  the  choice  between  the  00C  approach,  which  is 
rigorous  but  unefficient  from  an  economical  point  of  view,  and 
the  COC  approach,  which  is  computationally  inexpensive  but  some¬ 
times  unreliable.  In  order  to  reap  advantage  from  both  zero  and 
first  order  approximation  concepts,  a  hybrid  optimality  criterion 
(HOC)  was  developed  in  Refs,  f  1,31  ] .  In  the  HOC  approach  only 
.i  small  number  of  stress  constraints  are  treated  using  the  virtual 
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load  procedure  and  the  others  are  transformed  into  side  cons¬ 
traints  using  the  stress  ratio  formula.  The  selection  of 
constraints  requiring  first  order  approximation  can  be  made 
in  advance  on  the  basis  of  the  physical  judgement  of  the  de¬ 
signer.  It  can  also  be  performed  automatically  depending  on 
adequate  selection  criteria. 


A  simple  selection  criterion  consists  in  retaining  only 
the  potentially  active  stress  constraints  as  candidate  for 
first  order  approximation.  Potentially  active  constraints  are 
defined  as  those  which  are  close  to  their  limiting  values  at 
the  current  design  point  : 


a 


il 


-  o  . 

l 


(3.28) 


A  second  selection  criterion  states  that  a  stress  constraint 
can  be  treated  by  FSD  with  sufficient  accuracy  provided 


(3.29) 


that  is,  if  the  contribution  to  the  stress  comes  principally 

from  the  i  member  itself.  This  condition  arises  from  the  fact 

that  the  first  order  approximation  (3.5)  reduces  to  the  zero 

order  cne  (3.26)  in  a  statically  determinate  structure.  In  case 

of  structural  redundancy  both  approximations  are  thus  close  to 

each  other  if  condition  (3.29)  is  satisfied.  He ome t r i c a  1 1 y  this 

occurs  when  the  relevant  stress  constraint  for  the  ith  design 

variable  is  represented  in  the  design  space  by  a  surface  that  is 

t  h 

roughly  parallel  to  the  i  base  plane. 


both  s<  lection  criteria  (3.28)  and  (3.29)  must  of  course 
ho  repeated  after  each  structural  reanalysis.  They  must  be 
employed  within  some  adequate  tolerances  specified  by  the  user. 
Depending  upon  the  severity  of  these  tolerances,  the  number  of 
first  order  approximated  stress  constraints  will  be  small  or 
large  and  the  HOC  approach  is  accordingly  close  the  COC  or  to 
the  COC  approaches. 
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A  two-dimensional  geometrical  interpretation  of  the  COC, 

COC  and  HOC  approaches  is  given  in  Fig.  4,  in  the  space  of 
the  reciprocal  design  variables.  The  exact  stress  constraints 

o(  =  and  =  a ^  are  indicated,  as  well  as  their  linear 

.  .  %  — 

approximations  o  ^  l  =  °1  an“  °2  L  =  °2  tan8ent  at  T]  and  T2 

and  their  zero  order  or  FSD  approximations  ct,  _  =  a,  and 

-v  .  >>F  1 

o„  _  =  o  „  .  From  an  analysis  point  x  ,  solution  of  the  linea¬ 
ge  ,  r  2  o 

rized  problem  yields  the  point  C,  which  after  restoration  on 
the  true  constraint  surface  comes  in  C'  (COC  approach).  The 
exact  solution  lies  in  R.  The  FSD  approximation  leads  to  point 
C  and  after  scaling  to  point  C'  (COC  approach).  Finally  em¬ 
ploying  the  first  order  approximation  for  o  ^  and  the  zero  order 
one  for  a generates  the  point  H  and  after  scaling  the  point 
H'  (HOC  approach). 

3.3.  Thermal  and  body  forces 


In  this  section  the  concepts  developed  in  Ref.  f  I  1  are 
extended  to  take  into  account  applied  loads  that  depend  linearly 
on  the  design  variables,  for  example,  thermal  loads,  gravity 
field  loads,  etc...  The  static  structural  response  is  then  go¬ 
verned  by  the  relations 

Kq  =  g  =  g£  +  rv  (3.30) 


a  n  d 


o  =  Tq  -  o 
M  o 


(3.31) 


In  these  expressions  it  is  assumed  that  the  loads  g  can  be 
split  in  two  contributions,  namely,  the  external  loads  gc, 
which  are  independent  of  the  design  variables  (the  subscript 
c  stands  for  "constant”)  ,  and  the  "variable"  loads  g^ ,  which 
are  given  by  : 

n 


(3.32) 


40 


where  g  .  denote  constant  vectors.  The  initial  stresses  a 
vi  o 

are  ass ume d  to  be  due  to  thermal  effects  acting  without  any 
other  external  or  internal  loads  (residual  stresses).  Note 
that  the  index  l  is  omitted  in  this  section,  however,  the  re¬ 
sults  presented  can  readily  be  extended  to  take  into  account 
multiple  load  conditions. 

3.3.1.  The  scaling  concept 


The  key  to  being  able  to  construct  consistent  first  and  zero 
order  explicit  approximations  of  the  behavior  constraints  again 

lies  in  the  concept  of  scaling.  If  a'  denotes  a  point  in  the 

2 

design  space,  the  scaling  by  a  factor  f  yields  a  new  point  a 

a?  =  f  a!  (3.33) 

l  l 


at  which  the  displacements  and  stresses  become 


(3.34) 


and 


(3.35) 


where 


are  such  that 


(3.36) 

(3.37) 


qc  + 


(3.38) 


The  foregoing  expressions  are  apparent  form  the  fact  that 
the  displacements  associated  with  the  variable  loads 


art*  not  effected  by  scaling,  while  the  di  sp  lacements  q  due  to 

the  external  loads  g  must  be  divided  by  the  scaling  factor 

c 

(see  Ref.  (  I  ]  ,  section  3.1).  The  same  arguments  hold  for  the 
stresses  ;  only  the  contribution  is  modified  when  scaling 
i s  performed . 


Considering  the  flexibility  constraints  (1.6, 3.1)  and  the 
stress  constraints  (1.4, 3. 4)  it  can  be  seen  that  the  scaling 
factor  necessary  to  bring  the  design  point  back  on  the  res¬ 
traint  surfaces  are  given  respectively  by  : 


and 


where  it  is  understood  that 


c  J 


u  .  -  u  . 
J  V  j 


ck 


0.  -  o  .  +  o  . 

k  vk  ok 


(3.39) 


(3.40) 
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(3.42) 


3.3.2.  First  order  approximation 

In  order  to  construct  first  order  explicit  approximation  of 
the  behavior  constraints,  it  is  necessary  to  determine  analyti¬ 
cally  their  gradients.  For  that  purpose  the  pseudo-loads  tech¬ 
nique  discussed  in  section  3.1.2  will  be  employed.  Pifferen- 
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tiating  the  equilibrium  equations  (3.30)  furnishes  the  first 
partial  derivatives  of  the  nodal  displacements  as  the  solution 
of  the  linear  system  : 

K  Hr  ’  Hr  -  Irr  <  °-43) 

ill 


Since  the  stiffness  matrix  and  the  load  vector  are  both  linear 
in  the  design  variables  (see  Eqs .  1.9  and  3.32),  it  follows  that 


3a .  a  . 

l  l 


(3.44) 


3g  _  vi 


(3.45) 


so  that  by  solving  Eq.  (3.43)  : 


fa7  "IT  K'l<Ki  -  -  *vi> 

1  1 


(3.46) 


Therefore  the  gradients  of  the  flexibility  constraints  embodied 
in  K q  .  (1.6)  are  given  by 


i  1  T  - 1 

— =  -  —  b  !  K  (K  .  q  -  g  .  ) 

3a.  a .  j  i  vi 

l  l 

They  can  be  obtained  by  introducing  the  pseudo-load  vectors 


(3.47) 


*i  -  -  t:  (Ki<  ~  gvi)  -  -  it  <gi  "  gvi>  (3-48) 

i  i 

in  the  structural  analysis  phase.  Note  that  the  g  ^  '  s  represent 
the  loads  associated  with  the  i  element. 


At  this  point  it  is  important  to  notice  that  the  virtual 
load  technique  can  still  be  employed  to  evaluate  the  gradients 
of  the  flexibilities.  Using  again  virtual  load  cases  given 
numerically  hy  the  vectors  b ^  ,  it  appears  from  Eq  .  (3.47)  that 
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(3.49) 


with 


T  T 

c..  =  (q  K.  q.  -  g~.  q.)  a. 
lj  l  j  vi  1 J  l 


qj  (Ri  “  gvi)  ai  (3'50) 


where  q^  denotes  the  virtual  displacement  vector  just  as  in 
section  3.1.1.  It  is  evident  from  Eq  .  (3.49)  that  the  coef¬ 

ficients  c..  continue  to  be  equal  to  the  gradients  of  the 
flexibilities  with  respect  to  the  reciprocal  design  variables 
Therefore  the  first  order  explicit  approximation  of  the  fle¬ 
xibilities  can  now  be  constructed  by  using  first  order  Taylor 
series  expansion  : 
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which  reduces  to 

n 

n 

'V  o 

u  .  =  u  .  - 
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(3.52) 
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In  view  of  Eq.  (3.50) 

it  appears  that,  in  a  general  way 

: 
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ion  of  the 

flexibility  u. 

1  j 

has  been 

used 

is  this  equation  : 
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(3.54) 

u  .  = 

b.  c 

i  *  q  i 

K  q  =  q  K  q  . 

J 

J 

j 

J 

Finally,  the  first  order  explicit  approximation  of  a 
flexibility  constraint  turns  out  to  be 


u.  =u  .  +  c..  x. 

J  v .]  1=1  1  .1  1 


(3.55) 


where  the  uppercript  is  now  omitted  and  u^j  is  defined  as  the 
virtual  work  of  the  variable  loads  g^  on  the  virtual  displacements 

qi  * 
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(3.56) 


Or  course  the  same  developments  can  be  pursued  for  the  stress 
constraints,  yielding  the  first  order  explicit  approximations 


a  =  a  -a  +  .1  d.  x. 
k  vk  ok  i=l  lk  l 


(3.57) 


with 
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(3.58) 


and 


dik  "  (<1  Kiqk 


0;  ■  -  Evi)  •'*; 


’vi  k 


where  the  virtual  displacement  vector  is  due  to  a  virtual 
load  case  given  by  t  . 


With  regard  to  the  selection  of  either  the  pseudo-loads 
or  the  virtual  load  technique  for  computing  the  constraint 
gradients,  the  arguments  provided  in  section  3.1.3  are  still 
valid.  It  should  be  emphasized  that  again  the  virtual  load 
procedure  helps  understanding  the  nature  of  the  approximation 
made  in  linearizing  the  constraints.  Indeed  it  is  apparent 
from  Fqs.  (3.55  and  3.57)  that  the  points  of  intersection  of 
the  scaling  line  with  the  approximate  restraint  surfaces 

%  —  r\j  —  4  #  , 

u.  =  u.  and  a ^  are  precisely  given  by  the  scaling  factors 

f  and  f  obtained  in  Eqs.  (3.39  and  3.40)  for  the  real 

Uj  ok 

restraint  surfaces.  This  showns  that  Eqs.  (3.55,3.57)  continue 
to  represent  in  the  reciprocal  design  variable  space  the  tan¬ 
gent  planes  to  the  real  restraint  surfaces  at  their  points  of 
intersection  with  the  scaling  line.  In  fact,  the  quantities 

c .  . ,  u  d . .  and  (a  ,  -  a  ,  )  defining  the  first  order  explicit 

1.1  vj  lk  vk  ok 

approximations  of  the  flexibility  and  stress  constraints  (see 
Eqs.  1.55  and  3.57),  are  constant  in  the  case  of  a  statically 
determinate  structure  (the  explicit  problem  statement  is  then 


exact).  For  a  statically  indeterminate  structure,  they  depend 
implicitly  on  the  design  variables,  however,  they  remain  cons¬ 
tant  along  the  scaling  line. 


3.3.3.  Zero  order  approxima t i o n 

As  previously  defined, zero  order  explicit  approximation 
of  the  stress  constraints  is  obtained  by  applying  the  FSD 
criterion  (see  section  3.2).  It  corresponds  to  replacing 
the  stress  constraints  with  simple  side  constraints  resul¬ 
ting  form  the  stress  ratio  formula  (see  Eqs.  3.26,3.27). 

In  the  present  situation,  where  part  of  the  loads  dependsli- 
nearly  on  the  design  variables,  it  is  recommended  to  derive 
the  zero  order  approximation  from  the  first  order  one,  hy 
assumi ng  that 


where  6.,  denotes  the  kronecker  delta.  Eq .  (3.57)  then  reduces 
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(3.61) 


where  the  quantities  with  the  upperscript  °  are  frozen  to 
their  values  at  the  analysis  point.  The  stress  constraints 
are  therefore  transformed  into  the  side  constraints  (3.26) 
wi  th 


(3.62) 


With  this  new  definition  of  the  stress  ratio  formula,  the  ap- 

.  .  ,v 

proximate  constraint  surface  a  ^  =  a  ^  is  still  represented  by 

a  plane  normal  to  the  i axis  in  the  design  space  and  passing 

through  the  point  of  intersection  of  the  scaling  line  with  the 

true  constraint  surface  a.  =  a.. 

i  l 
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It  should  be  clearly  recognized  that  the  usual  stress 
ratio  formula  expressed  in  Eq .  (3.27)  for  multiple  load  cases, 

take  the  form 


(3.63) 


This  formula  is  equivalent  to  Eq .  (3.62)  only  in  the  case  of 

a  statically  determinate  structure  subject  to  mechanical 
and  thermal  loads,  because  of  the  following  identity  : 


a  .  -  a  .  =  0 
vi  oi 


(3.64) 


From  a  practical  point  of  view,  the  conventional  stress  ratio 
algorithm  (3.63)  is  simpler  to  employ  than  the  modified  one 
(3.62),  which  demands  separate  evaluation  of  the  pure  mecha¬ 
nical  stresses  o  .  and  the  other  contributions  (o  .  -  a  .). 

ci  o 1  v 1 
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4 .  DYNAMIC  CONSTRAINTS  ON  NATURAL  FREQUENCIES 

4.1.  Problem  statement 

The  problem  considered  in  this  chapter  consists  in  minimizing 
the  structural  weight  while  prescribing  lower  and  upper  limits 

on  the  natural  frequencies.  In  addition  side  constraints  are  as¬ 
signed  to  the  design  variables,  which  are  still  taken  as  the  trans¬ 
verse  sizes  of  the  structural  members.  Therefore  the  mathematical 
programming  problem  to  be  solved  reads  as  follows  : 


minimi ze 

w  = 

n 

P  . 

l  .  a  . 

(4.1) 

i  =  1 

l 

l  i 

subject 

to 
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l 

l 

The  dynamic 

cons  traints 

expressed  in 

Eq.  (4 

.  2 )  are 

directly  written 

in  te  rms  o  f 

the 

squares 

of 

the  frequencies. 

because 

these  quantities 

naturally  appear 

i  n 

the 

e i genp  r ob 1 em 

character i zing 

the  structural 

modal  analysis  : 


u)7 
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(4.4) 


In  this  equation,  K  and  M  respectively  represent  the  stiffness  and 
mass  matrices  and  (q.,  j  =  1,  m)  are  the  modal  displacements,  i.e., 

the  eigenvectors  solution  of  Eq .  (4.4  ),  associated  with  eigenvalues 

2 

0)  .  . 

J 

Several  approaches  have  been  proposed  to  solve  the  structural 
optimization  problem  embodied  in  Eqs.  (4.1  -  4.3).  Analytical  me¬ 
thods  were  employed  in  some  specific,  simple  situations  where  the 
problem  can  be  stated  in  the  form  of  differential  equations. 

NIORDSON  [  33  1  considers  transverse  vibrations  of  a  beam  simply 
supported  at  both  ends.  Using  variation  calculus  with  a  lagrangian 
multiplier,  he  maximizes  the  fundamental  frequency  for  a  given 
structural  weight.  Recoursing  to  a  similar  analytical  method, 
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TURNER  [34]  studies  axial  vibrations  of  a  rod  supporting  a  fixed 
mass.  The  problem  is  formulated  as  the  minimization  of  the  weight 
for  a  given  frequency.  TAYLOR  [  35  ]  reconsidered  the  same  problem 
by  employing  an  "energy"  approach  similar  to  that  of  PRAGER  [  36  ] 
for  static  constraints.  He  shows  that  the  problem  of  maximizing 
the  frequency  for  a  given  weight  yields  a  solution  that  also  renders 
minimum  the  weight  of  the  rod  for  a  given  frequency.  In  a  subse¬ 
quent  work,  TAYLOR  [  37  ]  introduced  a  lower  limit  on  the  rod  cross- 
sectional  area.  PRAOER  and  TAYLOR  [  38]  give  necessary  and  suffi¬ 
cient  optimality  conditions  for  achieving  maximum  frequency  at  given 
weight.  They  restrict  themselves  to  sandwich  structures  with  con¬ 
tinuously  varying  mass  and  stiffness  properties.  Using  a  similar 
energy  approach  SHEU  [39  ]  examines  one-dimensional  structures  such 
as  bars  and  beams,  in  the  numerically  interesting  case  where  mass  and 
stiffness  are  piecewise  constant. 

The  first  applications  of  mathematical  programming  techniques 
to  problems  involving  frequency  constraints  were  made  by  7.ARGHAMEE 
[  40  ] ,  who  employed  the  Rosen  gradient  projection  method  for  maximi¬ 
zing  the  fundamental  frequency  while  keeping  constant  the  structural 
weight  (linear  constraint).  Well  suited  to  complex  structural  sys¬ 
tems  analyzed  by  finite  elements,  this  design  procedure  needs  eva¬ 
luation  of  the  fundamental  frequency  gradient.  The  expressions  de¬ 
veloped  by  ZARGHAMEE  [  40  ]  are  exact  and  explicit  (see  next  sections). 
The  same  equations  were  also  derived  by  FOX  and  KAPPOOR  [32,41  1  and 
by  ROGERS  (  4  2  J  -  They  can  be  traced  back  to  the  works  of  WITTRICK 
I  43  1  and  in  some  sense,  to  the  work  of  FRAEIJS  de  VEUBEKE  [  44  ]  . 

RUBIN  [  45  ]  established  a  two  phase  procedure  consisting  in  maximi¬ 
zing  the  fundamental  frequency  at  constant  weight,  and  then  minimi¬ 
zing  the  structural  weight  at  constant  frequency.  This  second  phase 
is  accomplished  as  long  as  the  fundamental  frequency  remains  close 
enough  to  its  limiting  value.  Note  that  this  method  cannot  easily 
be  extended  to  the  case  where  multiple  frequency  constraints  are 
prescribed.  In  a  more  general  approach,  FOX  and  KAPPOOR  [  46  1  make 
use  of  the  Zoutendijk  feasible  direction  method  in  order  to  minimize 
the  structural  weight  subject  to  various  dynamic  constraints.  In 
addition  to  natural  frequency  constraints,  dynamic  response  cons- 
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traints  can  be  treated.  ROMSTAD  [  47  ]  uses  the  method  of  approxi¬ 
mation  programming  of  GRIFFITH  and  STEWART  [  48  ]  .  Although  restric¬ 
ted  to  a  single  constraint  in  the  applications  presented,  ROMSTAD's 
formulation  can  readily  be  extended  to  multiple  constraints  invol¬ 
ving  several  natural  frequencies.  HAUG,  PAN  and  STREETER  [  49  I  pro¬ 
pose  a  multi-purpose  method,  which  is  capable  of  dealing  with  se¬ 
veral  types  of  constraints,  including  frequency  constraints.  They 
employ  a  steepest  descent  method  in  conjunction  with  a  linearization 
process  for  the  governing  equations. 

The  concept  of  optimality  criterion  for  dynamic  problems  appears 
the  first  time  in  the  work  of  YOUNG  and  CHRISTIANSEN  [  50  ]  .  They 
design  a  space  truss  by  assuming  that  optimality  is  achieved  provided 
the  fully  stressed  design  criterion  is  satisfied  when  normalized  mo¬ 
dal  displacements  are  prescribed.  PRAGER  and  TAYLOR  [  38  ]  introduced 
an  analytical  optimality  criterion  for  continuous  sandwich  structures. 
No  attempt  is  made  of  deriving  redesign  formulas  for  discretized  sys¬ 
tems.  VENKAYYA,  KHOT  and  BERKE  [  51  ]  apply  their  optimality  criteria 
techniques,  initially  conceived  for  static  constraints,  to  the  case 
of  frequency  constraints.  In  the  examples  presented  the  recursive 
redesign  formulas  are  restricted  to  problems  involving  only  one  fre¬ 
quency  constraint.  VENKAYYA  and  KHOT  I  52]  subsequently  included 
cons  traints  on  dynamic  stresses  and  displacements.  TAIG  and  KERR  |*24],in 
an  innovative  approach  precursor  of  dual  methods,  considered  problems 
with  several  frequency  constraints.  Their  optimality  criteria  tech¬ 
nique  is  successfully  applied  to  realistic  structural  models  involving 
large  numbers  of  elements.  KAMAT  and  SIMITSES  [  53  ]  studied  beam 
models  lying  on  a  continuous  elastic  foundation  and  on  elastic  end 
supports.  They  propose  an  iterative  redesign  procedure  similar  to 
that  used  in  conventional  optimality  criteria  for  static  constraints. 
The  problem  statement  consists  in  maximizing  the  fundamental  frequen¬ 
cy  for  a  given  weight  of  the  beam.  Numerous  applications  are  pre¬ 
sented,  which  differ  from  each  other  according  to  the  support  condi¬ 
tions  and  the  dead  mass  distribution. 

In  the  present  work  the  basic  concepts  established  in  Ref.  I  I  1 
and  recalled  in  the  previous  chapters  will  be  followed,  that  is,  the 
original  problem  stated  in  Eqs.  (4.1  -  4.3)  will  be  transformed  into 
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a  sequence  of  explicit  problems.  The  explicit  problem  statement 
results  from  linearizing  the  frequency  constraints  with  respect 
to  the  reciprocal  design  variables.  Solving  partially  each  li¬ 
nearized  problem  using  a  primal  solution  scheme  yields  the  mixed 
method.  Solving  exactly  each  linearized  problem  using  a  dual  so¬ 
lution  scheme  leads  to  the  generalized  optimality  criteria  approach. 

4.2.  Application  of  the  mixed  method  concept 

As  in  the  case  of  static  stress  and  displacement  constraints 
considered  in  Ref.  [  1  ] ,  we  still  restrict  ourselves  to  structural 

models  made  up  of  bar  and  membrane  elements.  Therefore  the  stif¬ 
fness  and  mass  matrices  can  be  written 


K  -  i£. 


a  .  K . 

l  i 


(4.5) 


. I ,  a  .  M.  +  M 
1=1  1  1  c 


(4.6) 


where  K^,  and  are  independent  of  the  design  variables  a^ . 

and  denote  respectively  the  stiffness  and  mass  matrices  of 
the  if*'  element  when  a^  =  1  .  represents  the  contribution  of  the 

non- s t rue tur a  1  masses,  such  as  equipments,  fuel,  etc... 

Tn  order  to  implement  the  mixed  method,  consideration  must  first 
be  given  to  analytically  expressing  the  gradients  of  the  frequency 
constraints.  Differentiating  Eq.  (4.4)  with  respect  to  the  design 
variables  yields 


2  3M 


V  UU  •  A  U  • 

-  - i  M)  q.  +  (K  -  (0.  M)  - i 


(4.7) 


Pro multiplying  this  equation  by  and  taking  account  of  the  symmetry 

of  the  K  and  M  matrices  furnishes  : 
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(4.8) 


Since,  in  view  of  Eqs.  (4.5)  and  (4.6), 


(4.9) 


(4.  10) 


one  finally  obtains  the  gradients  under  the  form 

3uj.  <  I.*  M  ^ 

- — i  =  — -  q.  (K  .  -  a).  M.)  q. 

3ai  qTMq.  J  1  J  1  J 

J 


(4.11) 


In  order  to  facilitate  subsequent  developments,  let  us  introduce 
the  following  notations  : 


k..  =  (q..  K.  q . . ) a . 

ij  Mij  l  Mij  l 


(4.12) 


m..  =  (q..  M.  q..)a. 
i.l  ij  i  ij  i 


(4.13) 


y  .  =  q  .  M  q  . 

3  3  3 


n  m .  . 
13 


i.  =  q.M  q  .  =  y  .  -  .  I  .  — — 

3  3  c  3  i=l  a 


(4.14) 


(4.  15) 


where  q..  contains  the  components  of  the  j  eigenmode  associated 
13 

with  the  ic  element  ;  K.  and  M.  are  the  element  matrices,  whose 

1  1  t  ll 

dimension  equals  the  number  of  degrees  of  freedom  for  the  i  ele¬ 
ment.  The  gradients  obtained  in  Eq.  (4.11)  can  then  be  written  in 
the  compact  form 


„  2  2 

3w  .  k..  -  a),  m.. 

_ 3  =  13  3  13 

3a.  2 

l  y  .  a  . 
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(4.16) 


The  gradients  are  numerically  available  as  soon  as  the  modal  analysis 
has  bp en  performed,  yielding  the  natural  frequencies  u> j  and  the  cor¬ 


responding  eigenmodes  q. 


The  gradients  of  the  constraints  with  respect  to  the  reciprocal 
design  variables 

x  =  -L  (4.17) 


i  a . 
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which  are  of  primary  interest  in  this  study,  are  expressed  as 


-  2 
3t0  . 
J 


k  .  .  -  to  .  m .  . 


c  _ _ U  J  XJ 

3x.  ij  Pj 


(4.18) 
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The  explicit  forms  of  the  quantities  toj  can  now  be  obtained 


using  Taylor  series  expansion  around  the  design  point  x°  where 
the  modal  analysis  is  made  : 

2 

°  9 

(4. 19) 
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Restricting  the  expansion  to  the  first  order,  the  linearized 
eigenvalues  read  as  follows  : 


%2  2  11 

toT  =  to.  +  .  E ,  c?  .  (x.  -  x?) 
J  J  i=l  ij  i  i 


(4.20) 


where,  for  sake  of  convenience,  the  c..'s  denote  the  gradients 
of  the  eigenvalues  [see  Eq.  (4.18)1. 


Taking  account  of  the  definitions  (4.12)  and  (4.13),  it  can 
be  observed  that,  in  a  general  way 
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Finally  the  Taylor  series  expansion  given  in  Eq.  (4.20) 
becomes 
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From  the  foregoing  developments,  it  appears  that  the  linea¬ 
rized  forms  of  the  constraints  embodied  in  Eq.  (4.2)  are  ex¬ 
pressed  as  : 


n 

d!  <  . E  c?  .  x.  <  d° 
~J  i=l  iJ  i  J 


with 


2  2  m  i 
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(4.25) 
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Therefore  the  linearized  problem  upon  which  is  based  the  mixed 
method  can  be  generated,  by  replacing  the  real  constraints  (4.2) 
by  their  linearized  forms  (4.25)  : 

n  i  .  i  . 

minimize  W  =  .E  — - — —  (4.28) 

i=I  x . 
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subjectto 

d.  <  .1.  c.  .  x.  <  d.  j  =  1  ,  m  (4.29) 

-J  i=l  ij  i  j 


x.<  x.  <  x .  i  =  1  ,  n 

-l  l  l 


(4.30) 


Note  that  the  problem  is  defined  in  terms  of  the  reciprocal 
design  variables  and  that  the  index  0  is  now  omitted. 

The  linearized  problem  stated  in  Eqs.  (4.28-4.30)  exhibits 
the  same  characteristics  as  the  problem  L  constructed  in  chapter 
3  for  the  case  of  static  constraints  on  stresses  and  displacements. 
Therefore  the  concept  of  mixed  method  can  be  applied,  which  con¬ 
sists  in  partially  solving  problem  (4.28-4.30)  using  a  primal 
projection  method.  As  explained  in  section  2.2,  the  convergence 
control  parameter  k,  which  is  the  number  of  one -d i mens i ona 1  mi¬ 
nimizations  performed  at  each  stage,  must  be  chosen  according  to 
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the  degree  of  nonlinearity  that  the  frequency  constraints  exhibit 
in  the  reciprocal  design  space.  The  more  non  linear  the  cons¬ 
traints  are,  the  smaller  the  parameter  k  must  be  adopted.  It 
should  be  emphasized,  however,  that  even  in  the  limiting  case 
k  =  1,  the  mixed  method  does  no  longer  behave  as  a  pure  primal 
minimization  method  applied  to  the  real  problem.  Indeed  the  res¬ 
toration  phase  using  scaling  of  the  design  variables  does  not  lead 
necessarily  to  a  feasible  boundary  point  as  it  is  the  case  when 
stress  and  displacement  constraints  are  considered. 


The  effect  of  scaling  depends  now  upon  mass  distribution  in 
the  structure.  According  to  the  amount  of  non  structural  mass, 
scaling  indeed  produces  different  results.  The  effect  of  scaling 
in  best  investigated  by  examining  the  Rayleigh  quotient 


u 


2 

i 


(4.31) 


and  considering  the  two  following  special  cases. 


(1)  There  is  no  fixed  masses  :  the  structural  mass  constitutes 
then  the  only  contribution  to  the  mass  matrix,  i.e.,  1*^  =  0 
in  Eq.  (4.6).  Since  the  mass  and  stiffness  matrices  are 
both  linear  and  homogeneous  in  the  design  variables  I  see 
Eqs.  (4.5)  and  ( 4 . 6  )J,  it  is  clear  that  scaling  does  not 
modify  the  eigenvalues  nor  the  associated  eigenmodes. 

Note  also  that  both  the  numerator  and  the  denominator  in 
the  Rayleigh  quotient  (4.31)  are  multiplied  by  the  scaling 
factor,  so  that  each  frequency  aij  remains  unchanged. 


(2)  The  contribution  of  the  non  structural  mass  is  predomi¬ 
nating  and  therefore  the  structural  mass  can  be  neglected, 

i.e.  Eq.  (4.6)  reduces  to  M  =  M  .  Looking  at  the  Rayleigh 

c  2  . 

quotient  (4.31)  it  can  be  seen  that  the  eigenvalues  in¬ 
crease  in  proportion  to  the  scaling  factor.  The  associated 
eigenmodes  are  also  modified  and  must  be  reevaluated. 
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In  the  intermediate  case  where  structural  and  non  structural 
masses  contribute  with  the  same  order  of  magnitude  to  the  mass 
matrix,  the  scaling  process  requires  a  new  complete  finite  ele¬ 
ment  analysis  for  determining  the  modified  frequencies  and  eigen- 
modes  . 


In  this  connection  it  should  be  noted  that,  in  view  of  the 

numerical  results  obtained  up  to  now,  important  non  structural 

masses  seem  to  have  a  beneficial  effect  on  convergence  stability 

in  the  overall  optimization  process.  This  can  be  attributed  to 

the  fact  that  in  the  limiting  case  (2)  where  non  structural  masses 
.  2 

are  dominant,  the  terms  m^  appearing  in  the  constraint  gra¬ 
dients  (Eq  .  4.18)  are  small  when  compared  to  the  terms  k.^  (weak 
kinetic  energy).  The  structural  optimization  problem  then  be¬ 
haves  just  as  in  the  static  case,  where  the  mixed  method  concept 
has  proven  to  be  fully  valid.  In  the  other  cases,  recourse  to 
reciprocal  variables  does  no  longer  reduce  necessarily  the  non¬ 
linear  character  of  the  constraints,  which  is  the  fundamental 
idea  in  the  mixed  method.  In  such  a  case  it  might  be  useful  to 
adopt  more  explicit  approximations  for  the  constraints,  such  as 
the  second  order  Taylor  series  expansions  suggested  by  MIURA  and 
SCHMIT  [54  ] .  Note  also  that  when  there  is  no  fixed  mass  [  case 
(1)] ,  the  problem  is  ill-conditioned,  because  the  structural  weight 
can  be  rendered  arbitrarily  small  without  changing  the  frequency 
spec  t rum . 


Finally  it  is  worth  mentioning  that  in  the  SAMCEF  program 
[15  ]  ,  another  restoration  procedure  is  provided  to  replace  that 
of  scaling  when  necessary.  The  new  restoration  phase  consists 
in  maximizing  the  most  violated  frequency  constraint  while  keeping 
the  structural  weight  constant.  This  alternative  option  is  called 
whenever  the  frequency  constraints  reveal  violated  within  some 
given  tolerance  after  completing  the  minimization  phase.  The  mi¬ 
nimization  phase  still  proceeds  by  partially  solving  the  linearized 
problem  pcsed  in  F.qs.  (4.28  -  4.30).  This  process  of  alternating 
weight  minimization  and  frequency  maximization  is  based  on  the 
method  of  PUBIN  [  45  ]  discussed  in  section  4.1. 


4.3.  Generalized  optimality  criterion 


The  generalized  optimality  criterion  results  from  writing 
the  KUHN-TUCKER  conditions  for  the  linearized  problem  embodied 
in  Eqs.  (4.28  -  4.30).  It  is  apparent  that  this  explicit  op¬ 
timality  criterion  takes  exactly  the  same  form  as  in  the  case 
of  static  constraints,  that  is,  it  is  given  by  Eqs  .  (  2 . 2  1 -»2 . 2  5  ) 

A  physical  interpretation  of  the  generalized  optimality 
criterion  is  provided  by  introducing  respectively  the  elastic 
energy 
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and  the  kinetic  energy 
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in  the  i  finite  element  for  the  j  vibration  mode.  In  view 
of  the  definitions  (4.12)  and  (4.13),  these  energies  are  given 
by  : 
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and 
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The  difference  between  these  two  energies  character  i  zes  the  op¬ 
timality  of  each  structural  member.  Indeed  the  optimality  cri¬ 
terion  stated  in  Eq.  ( 2 . 2  1 )  can  also  be  written 
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where  the  rj's  are  the  lagrangian  multipliers  associated  with 
the  frequency  constraints. 
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The  criterion  reads  thus  just  as  in  the  static  case  (see 
section  3 . 1  .  l ) provided  that  the  concept  of  virtual  strain 
energy  is  replaced  with  that  of  difference  between  elastic 
energy  and  kinetic  energy.  In  the  special  case  where  only 
one  constraint  is  assigned  -  e.g.  lower  bound  on  the  funda¬ 
mental  natural  frequency  -  the  optimality  criterion  states 
that  in  the  optimal  structure  the  difference  between  the  elas¬ 
tic  energy  density  and  the  kinetic  energy  density  is  the  same 
in  each  element.  Note  that  this  condition  only  applies  to  the 
so-called  active  design  variables  (see  Eqs . 2 . 2  1 -2 . 23)  .  This 
optimality  criterion  was  first  proposed  by  PRAGER  and  TAYLOR 
I  38  1  in  the  context  of  continuous  models  (sandwich  structures). 
It  was  also  found  out  by  several  reseachers,  using  various  ap¬ 
proaches  [  34,40,45,511  • 


With  regard  to  the  redesign  formulas  the  general  procedure 
previously  established  in  the  static  case  remains  entirely  ap¬ 
plicable.  After  each  structural  analysis  phase,  yielding  the 
frequencies  and  the  corresponding  eigenmodes  q^  ,  the  coef¬ 
ficients  c„  are  computed  from  Eq .  (4.18).  The  explicit  problem 

posed  by  Eqs.  (4.28  -  4.30)  is  then  formed  and  it  is  solved  using 
either  primal  or  dual  algorithms  (see  sections  2.3,  2.4).  It 
should  be  clearly  recognized  that  solving  the  explicit  problem 
(4.27  -  4.29)  is  equivalent  to  finding  the  design  variables  that 
satisfy  the  optimality  conditions  (2.21-2.25) . Af ter  this  optimiza¬ 
tion  phase  is  completed,  a  structural  reanalysis  is  performed 
with  the  new  design  variables  and  the  process  is  repeated  until 
convergence  is  achieved.  It  is  important  to  point  out  that  in 
the  context  of  an  optimality  criteria  approach,  dual  algorithms 
are  better  than  primal  algorithms,  because  most  of  the  time  only 
a  small  number  of  frequency  constraints  must  be  taken  into  con¬ 
sideration.  On  the  other  hand,  if  it  is  found  necessary  to  con¬ 
trol  the  convergence  of  the  overall  optimization  process,  primal 
algorithms  are  obviously  the  most  appropriate  (mixed  method). 
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5  .  STABILITY  CONSTRAINTS  ON  LINEAR  BUCKLING  LOADS 

5.1.  Problem  statement 

The  structural  optimization  problem  for  the  case  of  global 
stability  constraints  on  the  linear  buckling  load  factors  exhi¬ 
bits  a  strong  analogy  with  the  problem  discussed  in  chapter  4, 
regarding  frequency  constraints.  The  mathematical  programming 
statement  of  that  problem  reads  indeed  as  follows 

n 

minimize  W  =  .1,  P  .  1  .  a  .  (5.1) 
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Just  as  the  natural  frequencies,  the  critical  load  factors 
are  defined  through  an  eigenproblem  : 


X  . 
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K  q  .  -  X  .  K  q  . 
J  J  G 


(5.4) 


where  K  represents  the  usual  stiffness  matrix  and  K  ,  the  geo- 
metrical  stiffness  matrix  ;  (q^  ,  j  =  I ,m)  denote  the  eigen¬ 

vectors  solution  of  problem  (5.4),  associated  with  eigenvalues 

X..  The  physical  meaning  of  the  q.'s  is  that  of  displacements 
J  J 

in  the  j  buckling  mode,  for  a  critical  load  factor  X^.  In 

the  rest  of  this  chapter  X.  will  be  simply  referred  to  as  "bu¬ 
ckling  load". 


The  first  attempts  of  optimizing  structures  subject  to 
stability  constraints  were  made  in  the  context  of  column  pro¬ 
blems,  with  the  aim  of  maximizing  the  buckling  load  for  a  gi¬ 
ven  structural  weight.  The  pioneering  work  of  CLAUSEN  [  55  ] 
was  concerned  with  the  optimal  shape  of  a  simply  supported  co¬ 
lumn  ;  more  recently  KELLER  [55  ]  reconsidered  the  same  problem. 

T APJBAKHSH  and  KELLER  [  57  1  extended  the  results  of  KELLER  [56  1 


to  the  case  of  columns  with  other  support  conditions  ;  they 
give  a  proof  of  optimality  for  their  solution.  In  their  ge¬ 
neral  structural  optimization  theory,  based  on  energy  consi¬ 
derations,  TAYLOR  [  58  ] ,  and  subsequently  PRAGER  and  TAYLOR 
[  38  ]  formulated  optimality  conditions  for  problems  invol¬ 
ving  maximization  of  the  critical  buckling  load.  TAYLOR  [  58  ] 
is  solely  concerned  with  columns  for  which  the  moment  of  i- 
nertia  is  proportional  to  the  square  of  the  cross-sectional 
area  ,  while  PRAGER  and  TAYLOR  [  38  ]  adopt  a  linear  relation 
between  moment  of  inertia  and  cross-sectional  area  (sandwich 
structures).  In  this  latter  work,  as  well  as  in  the  work  of 
TAYLOR  and  LIU  f  59  ] ,  minimal  allowable  values  are  assigned 
to  the  design  variables.  Finally  the  influence  of  self-weight 
on  the  optimum  design  was  examined  by  KELLER  and  NIORDSON  [  60  ] 
HUANG  and  SHEU  [61  ]  also  attacked  the  problem  of  self-weight 
in  optimal  columns. 

All  the  previously  mentioned  approaches  are  based  upon 
consideration  of  a  one-dimensional  continuous  model  and  they 
reduce  to  finding  a  single  variable  function.  The  discretized 
models,  to  which  the  finite  element  method  is  applicable,  were 
studied  by  ZARGHAMEE  [  62  1  »  in  the  context  of  a  mathematical 
programming  approach.  The  aim  is  to  minimize  the  structural 
weight  while  imposing  a  lower  bound  on  the  fundamental  buckling 
load.  To  this  end  ZARGHAMEE  [  62  1  employs  a  linearization  me¬ 
thod  related  to  the  cutting  plane  technique. 

The  optimality  criteria  approaches  are  well  suited  to  pro¬ 
blems  governed  by  buckling  phenomena,  because  in  addition  to 
the  usual  side  constraints,  only  one  behavior  constraint  has 
usually  to  be  taken  into  account  (lowest  critical  load).  In 
this  connection,  SIMITSES,  KAMAT  and  SMITH  [  63  1  propose  re¬ 
cursive  redesign  formulas  based  on  optimality  criteria  for  the 
case  of  elastically  supported  columns  subjected  to  distributed 
axial  loading.  The  finite  element  method  is  employed  for  ana¬ 
lyzing  the  discretized  model.  The  redesign  relations  are  si¬ 
milar  to  those  derived  by  BERKE  [  28  1  for  a  single  displacement 


constraint.  KIUSALAAS  [  64  ]  considers  a  more  general  problem 
where  two  buckling  modes  can  be  simultaneously  critical.  He 
employs  conventional  redesign  relations  for  multiple  constraints 
but  he  introduces  a  relaxation  factor  in  order  to  allow  for  con¬ 
trolling  the  convergence  of  the  optimization  process.  Results 
are  concerned  with  columns  and  portal  frames.  KHOT,  VENKAYYA 
and  BERKE  [  65  ]  »  following  out  the  general  optimality  criteria 
theory  previously  proposed  in  Ref.  [  51  ]  ,  developed  a  redesign 
procedure  for  minimizing  the  weight  of  complex  structural  systems 
subjected  to  a  single  buckling  constraint.  The  method  is  suc¬ 
cessfully  applied  to  portal  frames  and  trusses  involving  up  to 
132  finite  elements.  In  this  work  special  techniques  for  con¬ 
trolling  the  convergence  were  investigated,  such  as  over-rela¬ 
xation  and  generation  of  intermediate  design  points.  This  study 
was  motivated  by  instability  in  convergence  encountered  in  some 
problems  with  buckling  constraints. 

5.2.  Application  of  the  mixed  method  concept 

In  this  section  the  mixed  method  statement  for  problems  in¬ 
volving  buckling  constraints  will  be  based  upon  linearization 
of  the  constraints  in  the  reciprocal  design  variable  space,  just 
as  in  the  case  of  static  and  dynamic  constraints  previously  des¬ 
cribed.  This  requires  computation  of  the  critical  load  gradients, 
which  in  turn  assumes  to  be  known  the  geometrical  stiffness  ma¬ 
trix  dependence  on  the  design  variables.  For  the  structural  mo¬ 
dels  considered  herein,  the  stiffness  matrix  K  is  linear  in  the 
design  variables,  as  indicated  by  Eq.  (4.5).  This  is  not  the 
case  of  the  geometrical  stiffness  matrix  K  ,  for  which  the  expli- 
cit  dependence  with  respect  to  the  design  variables  disappears. 
The  matrix  K  can  still  be  split  into  the  contributions  of  the 
structural  members 
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where  K^,  denotes  the  geometrical  stiffness  matrix  of  the  i 
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(5.5) 


6  1 
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element.  The  element  matrix  K  is  related  to  the  initial  stress 

G  . 

state  in  that  element  and  therefore  it  depends  implicitly  on  all 
the  design  variables.  It  should  be  emphasized  that  the  geometri¬ 
cal  stiffness  matrix  is  independent  of  the  design  variables  for  a 
statically  determinate  structure. 


Differentiating  Eq .  (5.4)  with  respect  to  the  design  variables 

yields  : 
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T  .... 

Premul tipi  ted  by  q.  this  relation  furnishes,  just  as  in  the 

dynamic  case  (section  4.1)  : 
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From  Eq.  (5.5)  it  follows  that 
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In  opposition  with  the  static  and  dynamic  cases  previously  dis¬ 
cussed,  the  derivatives  appearing  in  Eq .  (5.8)  are  not  directly 

available,  because  the  elements  of  the  geometrical  stiffness  ma¬ 
trices  are  linear  functions  of  the  stresses  acting  in  the  pre¬ 
buckling  state.  Therefore  the  evaluation  of  the  derivatives  (5.8) 
requires  that  of  the  quantities 


(5.9) 


where  T  denotes  the  stress  matrix  of  the  kC^  element  and  q, 
are  the  generalized  displacements  of  the  k  element  resulting 
form  the  static  analysis  that  gives  the  prebuckling  state.  The 
stress  derivatives  (5.9)  can  be  computed  according  to  the  virtual 
load  or  pseudo-loads  techniques  discussed  in  section  3.1.  Note 
that  these  techniques  always  demand  introduction  of  a  cer^gin  num¬ 
ber  of  additional  fictitious  loading  cases.  Each  matrix  i R 
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then  constructed  in  the  same  way  as  the  corresponding  matrix 

K_,  ,  by  simply  replacing  the  stresses  o,  by  their  derivatives 
G.  K 

(5*9) . 


This  constitutes  of  course  a  costly  procedure,  which  for- 
tunetely  can  be  avoided  provided  the  terms  y—  appearing  in  the 
gradients  expressions  (5.7)  are  negligible.  This  assumption, 
which  is  typical  of  optimality  criteria  approaches  for  static 
constraints,  amounts  to  not  taking  into  account  the  effects  of 
structural  redundancy.  Adopting  this  assumption  it  follows  that 
the  gradients  of  the  critical  loads  can  be  written  in  the  form  : 
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(5.10) 


where  the  coefficients  k..  are  given  by 
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and 


(5.12) 


From  a  practical  point  of  view,  the  gradients  are  numerically 
available  as  soon  as  the  eigenproblem  (5.4  )  has  been  solved, 
yielding  the  critical  loads  X^  and  the  associated  buckling  mo¬ 
des  q  j .  A  preliminary  static  analysis  must  be  accomplished  to 
permit  evaluation  of  the  geometrical  stiffness  matrix  K^,. 


When  the  now  well  known  change  of  variables  (4.17)  is  intro¬ 
duced, the  gradients  of  the  critical  loads  become 
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where  the  notation  c^j  is  again  employed  for  sake  of  conciseness. 
By  pursuing  the  usual  linearization  process  (see  sections  2. 1,3.3 
and  4.2),  it  come s 
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where  x°  denotes  again  the  reciprocal  design  point  at  which 
structural  analysis  is  performed.  Taking  into  account  the 
definitions  (5.11)  and  (5.13),  it  can  be  observed  that 
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As  a  consequence  the  linearized  forms  of  the  critical  loads 
developed  in  Eq.  (5.14)  reduce  to 

^  .  =  2  X?  +  .  E  c?  .  x.  (5.16) 

J  J  ij  i 

The  buckling  constraints  stated  in  Eq.  (5.2),  when  linearized 
at  the  design  point  x° ,  read  then  as  follows  : 

n 

X.  -  2  X?  <  .1.  c° .  x.  <  X.  -  2  X?  (5.17) 
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The  notion  of  mixed  method  is  directly  applicable  to  the 

problem  considered  in  this  section.  The  scaling  of  the  design 

variables  constitutes  indeed  an  exact  restoration  procedure.  As 

in  the  static  case  (i.e.  with  static  stress  and  displacement 

constraints),  scaling  does  not  demand  any  structural  reanalysis. 

When  all  the  direct  design  variablesare  multiplied  by  the  same 

number,  the  critical  buckling  loads  are  also  multiplied  by  this 

scaling  factor.  However  the  scaling  operation  has  now  an  effect 

upon  the  value  of  the  coefficients  c^.  In  opposition  with  the 

static  case,  where  the  coefficients  c..  remain  constant  along  the 

i  J 

scaling  line,  they  are  here  multiplied  by  the  square  of  the  sca¬ 
ling  factor.  In  other  words,  if 
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then  it  follows  that 
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In  these  equations  the  upperscript  °  stands  for  the  initial  de¬ 
sign  where  structural  analysis  is  made  and  the  upperscript  1  re¬ 
fers  to  the  scaled  design,  f  denoting  the  scaling  factor. 


As  a  result,  in  the  space  of  the  reciprocal  design  variables, 
the  linearized  expressions  given  by  Eq.  (5.16)  furnish  the  exact 
values  of  the  critical  loads  and  of  their  gradients  only  at  the 
design  point  x°  where  the  structural  analysis  is  performed. 

These  linearized  forms  are  no  longer  correct  up  to  the  first 
order  at  any  point  along  the  scaling  line  passing  through  x°,as 
it  was  the  case  for  static  stresses  and  displacements.  Geome¬ 


trically  it  can  be  seen  that  the  real  restraint  surfaces  A ^  =  A^ 
and  Aj  =  A^  are  represented  in  the  reciprocal  design  variable 
space  by  the  planes  given  by  Eq .  (5.17)  (with  equality  signs). 


Such  a  plane  does  no  longer  pass  through  the  point  of  intersec¬ 


tion  of  the  corresponding  real  restraint  surface  with  the  sca¬ 


ling  line  (see  Fig.  1  ). 


However,  by  making  use  of  the  expressions  (5.19,  5.20)  it  is 

possible  to  relate  the  critical  loads  A^  and  their  gradientsc^j 

along  the  scaling  line,  without  having  to  recourse  to  additional 

structural  reanalyses.  Therefore  tangent  planes  to  the  restraint 

surfaces  can  still  be  constructed,  by  linearizing  each  buckling 

constraint  at  its  point  of  intersection  with  the  scaling  line, 

rather  than  directly  at  the  design  point  where  the  structural 

analysis  is  accomplished.  From  Eq  .  (5.19)  it  can  be  concluded 

that  the  scaling  factors  necessary  to  bring  the  design  point  back 

to  the  real  restraint  surfaces  A.  =  A.  and  A.  =  A.  are  respecti- 
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Consequently  in  view  of  the  linearized  forms  (5.16),  the  equations 

of  the  tangent  planes  to  the  surfaces  X  .  =  _X •  and  A.  =  X.  are  res- 
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pectively  given  by 
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From  the  foregoing  developments  it  is  easily  verified  that 
by  setting 
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the  linearized  problem  upon  which  is  based  the  mixed  method  concept, 
reads  as  follow  : 
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(5.28) 

Because  of  the  way  this  explicit  problem  has  been  generated,  the 
mixed  method  can  be  employed  exactly  like  in  the  static  case 
(section  2 .  3’>  .  However  it  should  be  clearly  recognized  that  the 
exolicit  expressions  embodied  in  Eq.  (5.14)  are  the  true  lineari¬ 
zed  forms  of  the  critical  loads  only  if  the  structure  is  stati¬ 
cally  determinate.  Otherwise  they  contain  an  error  that  is  rela¬ 
ted  to  the  degree  of  structural  redundancy.  When  this  error  is 
negligible,  or  when  the  coefficients  c^.  are  computed  using  the 
exact  gradients  given  in  Eq.  (5.7),  then  the  linearized  problem 
stated  in  Eqs.  (5.26-5.28)  can  be  geometrically  interpreted  just 
as  the  problem  L  developed  in  section  2.1  for  the  case  of  stress 
and  displacement  constraints.  Namely,  each  real  restraint  surface 
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(Aj  =  A.  or  Aj  =  Aj)  is  replaced  in  the  reciprocal  design  variable 
space  with  its  tangent  plane  at  its  point  of  intersection  with  the 
scaling  line  (see  Fig.  1  ). 

In  these  circumstances,  the  mixed  method  described  in  section 
2.3  can  be  effectively  used  to  control  the  convergence  of  the  op¬ 
timization  process.  In  the  special  limiting  case  where  only  one 
minimization  step  is  completed  at  each  stage  (k  =  1),  the  mixed 
method  still  reduces  to  a  pure  projection  method  of  non  linear 
programming.  Recourse  to  the  control  parameter  k  should  in  fact 
reveal  especially  useful  for  optimization  problems  involving  glo¬ 
bal  buckling  constraints.  Indeed  KHOT,  VENKAYYA  and  BERKE  [  65  ] 
report  some  difficulties  in  achieving  convergence  when  applying 
their  optimality  criteria  techniques  to  such  problems.  The  mixed 
method,  which  lies  between  optimality  criteria  and  mathematical 
programming  approaches,  seems  capable  of  overcoming  these  diffi¬ 
culties.  Unfor tunetely ,  in  the  current  version  of  the  SAMCEF 
program  [  15  ] ,  in  which  the  mixed  method  is  available,  the  con¬ 
vergence  control  capability  has  not  been  tested  extensively  enough 
for  buckling  constraints  so  as  to  provide  a  definitive  answer  with 
regard  to  its  efficiency. 


5 . 3  Generalized  optimality  criterion 


As  explained  in  section  2.3,  the  generalized  optimality  cri¬ 
teria  approach  consists  in  solving  exactly  at  each  stage  the  li¬ 
nearized  problem  stated  in  Eqs.  (5.26-5.28).  Alternatively  the 
optimality  criterion  can  be  written  explicitly  by  using  the  KUHN- 
TIJCKFR  conditions.  Just  as  in  the  static  and  dynamic  cases,  these 
conditions  take  the  form  expressed  in  Eq  s  .(2 . 2 1 -2 . 25)  with  the  pro¬ 
per  definition  of  the  coefficients  C£j*  Another  form  of  the  op¬ 
timality  criterion,  which  is  more  suitable  to  physical  interpre¬ 
tation,  reads  as  follows  : 
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represents  the  strain  energy  in  the  i  element  for  the  j  bu¬ 
ckling  mode.  Note  that  each  lagrangian  multiplier  r.  must  be 
positive  or  zero  according  to  whether  the  associated  buckling 
constraint  is  active  or  not. 

The  criterion  possesses  an  interesting  physical  interpreta¬ 
tion  in  term  of  strain  energy  densities,  in  a  similar  way  as  for 
the  static  case  (section  3.1.1).  In  the  special  but  frequent  si¬ 
tuation  where  only  one  buckling  constraint  is  assigned,  the  op¬ 
timality  criterion  states  that  in  the  optimal  structure,  the 
strain  energy  density  associated  with  the  critical  buckling  mode 
is  the  same  in  each  element.  Several  researchers! 38  »  63  »  64  »  65 
have  proposed  the  optimality  criterion  posed  by  Eq  .  (5.29),  but 

without  deriving  efficient  redesign  relations. 

As  previously  reported  in  this  work  on  several  occasions, 
the  only  rigorous  way  of  finding  the  design  variable  values  that 
satisfy  the  optimality  criteria  equations  is  to  solve  exactly 
the  linearized  problem  using  efficient,  special  purpose  mathema¬ 
tical  programming  algorithms.  So  the  redesign  procedure  is  as 
follows  :  after  each  structural  reanalysis,  yielding  the  coef¬ 
ficients  c..j  via  Eqs.  (5.11,  5.13  and  5.25),  the  new  estimates 
of  the  optimal  design  variables  are  generated  by  applying  pri¬ 
mal  or  dual  algorithms  to  the  explicit  problem  embodied  in  Eqs. 
(5.26-5.28).  As  in  the  dynamic  case  considered  in  chapter  4, 
dual  algorithms  are  most  often  recommended,  because  the  number 
of  critical  buckling  constraints  is  usually  small. 
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6.  EXTENSION  TO  FLEXURAL  MODELS 


The  structural  optimization  methods  reported  in  Ref.  (  I  ] 
and  in  the  previous  chapters  of  this  work  were  essentially 
developed  in  the  context  of  thin-walled  structures  for  which 
finite  element  idealization  by  means  of  bars  and  membranes 
leads  to  adequate  estimation  of  the  quantities  describing  the 
structural  response  (e.g.,  stresses,  displacements,  natural 
frequencies,  critical  buckling  loads).  The  stiffness  and  mass 
matrices  of  these  structural  models  are  linear  in  the  design 
variables,  which  so  far  were  the  transverse  sizes  of  the  elements. 
This  situation  of  course  facilitates  derivation  of  the  mixed 
method  and  generalized  optimality  criterion. 

Such  a  finite  element  representation  is  valid  only  when 
each  structural  member  can  be  assumed  to  be  in  a  plane  state  of 
stress.  So  at  the  element  level  the  applied  forces  must  act 
mainly  in  extension.  If  the  flexural  loading  cannot  be  neglec¬ 
ted,  then  it  is  necessary  to  introduce  other  types  of  finite 
elements  than  the  bars  and  membranes  considered  up  to  now.  The 
stiffness  matrix  of  a  bending  element  such  as  a  beam  or  a  plate 
is  no  longer  proportional  to  its  transverse  size  and  therefore 
the  optimization  methods  proposed  in  the  previous  sections  must 
be  modified. 

6.1.  Choice  of  the  design  variables 

The  optimum  design  of  structures  for  which  the  stiffness 
properties  do  not  depend  linearly  on  the  design  variables  was 
mainly  investigated  in  the  case  of  continuous  models.  Most 
often  the  optimization  problem  is  restricted  to  constraints 
on  the  natural  frequencies  and  on  the  buckling  loads  [66,67,68]  . 

In  this  chapter  attention  is  first  focused  on  discretized  models 
made  up  of  pure  beam  and  plate  elements.  Next  consideration 
will  be  given  to  the  more  general  case  where  flexion  and  exten¬ 
sion  forces  are  simultaneously  applied  at  the  element  level. 
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6.1.1.  Beam  elements 

The  way  to  deal  with  structural  optimization  problems  for 
a  beam  subject  to  pure  bending  depends  upon  the  relationship 
between  the  principal  moment  of  inertia  I  and  the  cross-sec¬ 
tional  area  a.  A  wide  variety  of  situations  is  taken  into 
account  by  means  of  the  following  relation  : 

I  =  c  a*5  (6.1) 

where  c  is  a  constant  that  depends  only  on  the  shape  of  the 
beam  cross-section  and  p  is  a  positive  number  [  53,  64  ) .  Most 
of  the  time  p  is  taken  as  an  integer  number,  equal  to  1,  2  or  3. 

The  case  p  =  1  corresponds  to  thin-walled  beams,  for  example, 
sandwich  beams,  pipes  with  fixed  diameter  but  variable  thickness, 
etc...  The  methods  described  in  sections  3,  4  and  5  remain  then 
fully  applicable,  since  the  stiffness  matrix  of  each  structural 
member  continues  to  be  proportional  to  the  design  variable,  which 
can  be  either  the  cross-sectional  area  or  the  moment  of  inertia. 

In  the  context  of  discretized  models,  VF.NKAYYA  [  23  ]  was  the  first 
to  study  this  kind  of  problem.  He  makes  the  observation  that 
when  the  ratio  of  the  moment  of  inertia  to  the  cross-sectional 
area  is  kept  constant,  the  scaling  concept  remains  entirely  valid. 
This  is  an  important  property  from  the  point  of  view  of  optimality 
criteria  approaches 

The  case  p  =  2  is  that  of  beams  with  uniformly  varying  cross- 
section;  the  shape  of  the  cross-section  is  kept  constant  while  its 
area  varies  along  the  beam  axis  (dilatation  or  contraction).  This 
type  of  beam  problem  was  analytically  investigated  by  KELLER  (  56  ] 
and  TAYLOR  [  58  ]  ,  for  a  single  constraint  on  the  critical  buckling 
load  (Euler  beam). 

Finally  the  case  p  =  3  is  concerned  with  beams  with  full  cross- 
section  whose  height  varies  wttile  ether  sizes  are  fixed  (case  similar  to 
that  of  plate  elements)  (53  ]. 
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In  order  to  illustrate  the  versatility  of  relation  (6.1), 
let  us  consider  a  beam  with  rectangular  cross-section.  It 
follows  that 


I 


(6.2) 


and 


a  =  b  h 


(6.3) 


where  b  is  the  width  of  the  cross-section  and  h,  its  height. 
Varying  the  width  b  and  keeping  constant  the  height  h  leads 
to 


1  =  12  3 


(6.4) 


which  corresponds  to  c  =  y-y  and  p  =  1  in  Eq  .  (6.1)  (linear 


re lation) . 
that 


If  now  the  height  h  is  the  variable  ,  it  follows 


I  =  — j  a3  (6.5) 

1  2b 

namely,  c  =  — - — y  and  p  =  3  in  Eq.  (6.1)  (cubic  relation). 

1  2  b  ^ 

Finally  keeping  the  ratio  h/b  constant  (constant  shape  cross- 
section)  yields 


I 


that  is,  a  quadratic  relation. 


(6.6) 


On  the  other  hand,  in  the  case  of  complex  cross-sectional 
shapes,  the  number  p  appearing  in  Eq.  (6.1)  can  be  taken  non 
integer,  by  recoursing  to  empirical  formulas  relating  the  mo¬ 
ment  of  inertia  to  the  cross-sectional  area  for  standard  gage 
sizes.  MOSES  and  ONADA  [69  1  ,  as  well  as  KAVLIE  and  MOE  [  70  ]  , 
employ  such  formulas  for  optimizing  elastic  grids.  Note  for 
example  the  following  relation,  which  is  deduced  from  English 
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norms  for  beam  construction  : 


I 


2.65 


'•°07  (T7480> 


(6.7) 


This  expre s s ion  leads  to  the  values  c  =  0.356  and  p  =  2.65  in 
the  general  relation  (6.1). 


Once  the  structural  model  is  characterized  via  an  explicit 
relation  of  the  form  given  in  Eq .  (6.1),  the  problem  of  choosing 
appropriate  design  variables  disappears.  The  natural  choice 
for  expressing  the  weight  of  a  beam  is  that  the  main  design  va¬ 
riable  must  be  the  cross-sectional  area.  However  the  stiffness 
properties  are  primarily  dependent  on  the  moment  of  inertia  and 
therefore  it  can  be  expected  that  the  explicit  approximations 
of  the  behavior  constraints  will  be  the  most  accurate  when  written 
in  terms  of  the  reciprocal  moments  of  inertia.  Since  both  quan¬ 
tities  -  cross-sectional  area  and  moment  of  inertia  -  are  connec¬ 
ted  to  each  other  through  Eq .  (6.1),  either  one  can  be  taken  as 
independent  design  variable.  For  a  beam  subjected  to  pure  ben¬ 
ding,  the  flexural  rigidity  is  proportional  to  the  moment  of  i- 
nertia  and  therefore  in  a  finite  element  context  the  structural 
stiffness  matrix  exhibits  the  following  explicit  form  in  terms 
of  the  cross-sectional  areas  : 

n  p  _ 

K  =  a .  K .  p  >  0  (6.8) 

i=  1  l  l 

where  each  matrix  is  independent  of  the  design  variables  a^. 

6.1.2.  Plate  elements 


With  regard  to  plate  elements  subjected  to  pure  bending,  two 
cases  must  be  distinguished.  The  first  case  is  that  of  sandwich 
plates  with  constant  core  thickness.  Thesheet  thicknesses  cons¬ 
titute  then  the  design  vari  ables.Theref  ore  the  stiffness  matrix 
continues  to  depend  linearly  upon  the  design  variables  and  the 
methods  previously  developed  in  this  work  remain  fully  applicable. 
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The  second  case  is  concerned  with  conventional  plates  with 
variable  thickness.  For  this  type  of  structural  model  the 
stiffness  is  proportional  to  the  cube  of  the  thickness  and 
hence  the  relation  (6.8)  must  be  chosen  with  p  =  3.  In  the 
context  of  continuous  models,  plate  optimum  design  was  inves¬ 
tigated  by  OLHOFF  [  71  J  and  ARMAND  [  72  ] ,  for  the  case  of  a 
single  constraint  on  the  fundamental  natural  frequency,  and  by 
FRAUENTHAL  [  73  ] ,  for  the  case  of  specified  critical  buckling 
load. 

6.1.3.  Flexion  -  extension  elements 

When  flexion  and  extension  loadings  act  simultaneously  with 
comparable  intensity  at  the  element  level,  the  definition  (6.8) 
of  the  stiffness  matrix  can  no  longer  characterize  the  structural 
model  with  sufficient  accuracy.  In  order  to  illustrate  this  fact, 
let  us  consider  again  the  beam  with  rectangular  cross-section 
previously  discussed  and  assume  that  only  the  cross-section  width 
can  be  varied.  If  the  beam  is  subjected  to  combined  flexion, 
some  terms  of  the  stiffness  matrix  require  p  •  I  in  Eq.  (6.1) 
(principal  moment  of  inertia  proportional  to  the  design  variable), 
while  others  demand  p  =  3  (secondary  moment  of  inertia  proportional 
to  the  cube  of  the  design  variable). 


In  these  circumstances,  it  is  necessary  to  recourse  to  a 
more  complete  definition  of  the  element  stiffness  matrices  in 
terms  of  the  design  variables.  As  suggested  by  KIUSALAAS  [  64,74]  , 
the  following  expression  can  be  chosen  : 


K 


(6.9) 


where  the  matrices  K^P^  are  independent  of  the  design  variables 
a..  Each  matrice  must  be  memorized  independently  for  each 
finite  element,  in  such  a  way  that  the  contributions  of  the  va¬ 
rious  terms  can  be  subsequently  separated. 
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The  expression  (6.9)  is  rather  general  and  it  permits 
great  variety  in  the  choice  of  the  element  sizes  that  can 
be  varied.  As  an  example  let  us  consider  a  beam  element 
with  hollow  circular  cross-section  (thick-walled  pipe). 

The  extensional  stiffness  is  included  in  Eq.  (6.9)  through 
a  term  proportional  to  the  design  variable,  which  is  taken 
as  the  cross-sectional  area,  while  the  dependence  of  the  fle¬ 
xural  stiffness  upon  the  design  variable  is  related  to  the 
way  the  beam  cross-sectional  sizes  can  be  varied.  So  if  the 
cross-section  is  kept  similar  to  itself  all  along  the  length 
of  the  tube,  then  it  follows  that 

K.  =  K(.1  }  a.  +  K^2)  a ?  (6 .10) 

11111 


If  now  the  external  diameter  varies,  while  the  internal  one 
is  fixed,  then  another  relation  must  be  adopted  : 


K. 

l 


a.  ♦  K<3)  a3 
1111 


(6.11) 


Based  on  similar  arguments,  it  can  be  seen  that  a  flat  shell 
element,  made  up  of  the  assembling  of  a  membrane  element  (ex¬ 
tension)  and  a  plate  element  (flexion),  leads  to  two  distinct 
contributions  to  the  stiffness  in  Eq.  (6.9)  :  one  is  linear 
in  the  thickness  and  the  other  one  is  cubic,  just  as  in  Eq .  (6.11). 


On  the  other  hand  the  mass  matrix  of  each  finite  element 
continues  to  exhibit  a  linear  form  in  terms  of  the  design  va¬ 
riables  : 


M.  =  l  aP  (6.12) 

l  p»o  1  1 

By  following  basically  the  same  procedures  as  for  the  case 
of  thin-walled  structures  (see  chapters  3,  4  and  5),  it  is  pos¬ 
sible  to  construct  explicit  optimality  criteria,  provided  the 
contributions  of  the  various  element  matrices  (and  for 

frequency  constraints)  are  decomposed  in  the  constraint  explicit 
approx ima  t ions . 


6.2.  Application  of  the  mixed  method  concept 


From  the  foregoing  developments  it  follows  that  a  general 
automatic  redesign  procedure  can  be  generated,  by  following 
out  the  now  customary  scheme  for  deriving  the  mixed  method  and 
the  generalized  optimality  criterion  (see  section  3.3,  4.2  and 
5.2).  In  order  to  outline  the  procedure  followed,  attention  is 
restricted  to  the  case  of  constraints  on  the  static  structural 
response,  i.e.,  upper  limits  on  stresses  and  displacements.  It 
is  assumed  that  the  structural  discretization  is  made  up  entirely 
of  elements  of  the  same  type.  Therefore  the  stiffness  matrix 
exhibits  the  form  given  in  Eq .  (6.8),  where  p  takes  on  the  same 

value  for  each  member.  To  fix  ideas  the  method  development  is 
limited  to  displacement  constraints,  keeping  in  mind  that  stress 
constraints  can  always  be  converted  into  displacement  constraints 
(first  order  approximation)  or  replaced  with  minimum  gage  cons¬ 
traints  (zero  order  approximation). 

Just  as  in  the  case  of  th in-wa lied  structures,  the  mixed  me¬ 
thod  concept  relies  on  a  change  of  variables  tending  to  reduce 
the  nonlinear  character  of  the  constraints  : 


x . 
i 


(6.13) 


The  next  step  is  to  linearize  the  constraints  in  terms  of  the 
new  variables  x^,  which  requires  gradient  evaluation  (see  for 

example  Eq .  2.11).  Differentiating  the  equilibrium  equation 

.  .  .  .  T 

K  q  =  g  and  using  the  definition  of  a  flexibility  Uj  =  b^  q 

leads  to  the  constraint  gradients  via  the  pseudo-loads  techni¬ 
que  (see  section  3.1.2)  : 


(6.14) 


Alternatively,  introducing  the  virtual 
associated  with  the  flexibilities,  the 
written  as  : 


oad  vectors  q.  *  K  '  b. 

J  J 

;radients  (6.14)  can  also  be 
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T  3K 
lj  3x. 


(6.15) 


This  latter  expression  corresponds  to  the  definition  of  the 
gradients  using  the  virtual  load  technique  (see  section  3.1.1). 

In  terms  of  the  new  variables  x^,  the  global  stiffness  ma¬ 
trix  given  in  Eq.  (6.8)  becomes 


K  =  .1  .  — 
1  =  1  x 


(6.16) 


Its  first  partial  derivatives  are  therefore  : 


(6.17) 


when  K.  is  the  stiffness  matrix  of  the  i th  element.  Finally, 

l  ... 

by  substituting  Eq.  (6.17)  into  Eq .  (6.15),  the  flexibility 

gradients  take  the  form  : 


•  T  V 

-  q  .  K  .  q 

x;  J  1 


(6.18) 


so  that  the  following  identity  holds 


n  3u  . 

u.  =  .£.  (-^-i)  x. 

j  1=1  dX.  1 


(6.19) 


By  employing  the  general  linearization  process  outlined  in  sec¬ 
tions  3.3,  4.2  and  5.2,  we  can  express  the  explicit  forms  of 
the  flexibilities  as  : 


n  n 

a.  a  o  o  o 

u.=u.-.E,c..x.+.£.c..x. 
j  j  i=l  lj  i  t*i  tj  i 

where  the  coefficients  c°j  denote  the  flexibility  gradients 
evaluated  at  the  design  point  x°  : 


(6.20) 


C?.  =  (-T-J-) 

1  j 


(6.21  ) 
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In  view  of  Eq .  (6.19),  it  is  apparent  that  the  first  and  second 

terms  in  Eq .  (6.20)  cancel  each  other,  so  that  the  linearized 
constraints  exhibit  the  same  form  as  in  the  case  of  bar  and  mem¬ 
brane  elements  previously  Considered  (omitting  now  the  uppers- 
cr l p  t  )  : 


u.  =  .1,  c..  x.  <  u  .  j  =  l,m 

j  i=l  ij  l  j 


(6.22) 


Therefore  the  linearized  problem  statement  is  similar  to 
that  stated  in  Eqs.  (2. 15-2. 17),  except  that  the  objective  func¬ 
tion  must  be  replaced  by  : 

n  P.  I . 

"  '  iE-,  -T77  <6-23) 

X  . 

1 

The  linearized  forms  of  the  constraints  embodied  in  Eq .  (6.22) 
are  still  exact  up  to  the  first  order  at  any  point  along  the 
scaling  line  in  the  space  of  the  x^'s  (generalized  reciprocal 
space),  or  in  the  space  of  the  a?.  So  the  scaling  process  keeps 
all  its  properties  (see  sections  3.1  and  3.2),  provided  scaling 
is  performed  on  the  quantities  a?,  rather  than  directly  on  the 
design  variables  a^.  As  a  result  scaling  can  still  be  employed 
to  bring  the  design  point  back  on  the  boundary  of  the  feasible 
region  (restoration  phase).  The  mixed  method  thus  behaves  exac¬ 
tly  as  in  the  case  of  structural  models  made  up  of  bar  and  mem¬ 
brane  elements  (see  section  2.3).  The  parameter  k  continues  to 
permit  control  over  convergence  in  the  optimization  process. 

In  the  limiting  case  k  =  1,  the  mixed  method  identifies  itself 
to  a  primal  mathematical  programming  algorithm. 


6.3.  The  generalized  optimality  criterion 


As  explained  in  section  2.4.2  the  explicit  form  of  the  ge¬ 
neralized  optimality  criterion  is  produced  by  expressing  the 
KUHN-TUCKER  conditions  for  the  linearized  problem.  When  recast 
in  terms  of  the  direct  design  variables  a^,  the  linearized  pro¬ 
blem  takes  the  form 
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minimi ze 


W  =  .  Z  ,  p .  1.  a . 

i=l  l  l  l 


(6.24) 


subject  to 


n  c  .  . 

.  E  — —  <  u  . 

1=1  ap  J 
l 


j  =  1  ,m 


(6.25) 


a  .  <  a  .  <  a  . 
-l  l  l 


i  =  1  ,  n 


(6.26) 


The  optimality  criteria  equations  are  very  similar  to  that  em¬ 
bodied  in  Eqs.  (2.21-2.25),  provided  that  care  is  taken  of  the 
exponent  p  appearing  in  Eq .  (6.25).  For  example  the  redesign 
relations  for  the  active  design  variables  must  read  as  follows  : 


a  .  =( 


c.  .  r  . ) 


i  P  j/i  J  i J  J 


(6.27) 


where  it  is  understood  that  the  dual  variables  r^  (i.e.  the 
lagrangian  multipliers)  must  satisfy  the  complementarity  con¬ 
ditions  (2.24,  2.25). 

To  obtain  a  physical  interpretation  of  the  optimality  criterion, 
we  introduce  the  virtual  strain  energies 


e  .  .  =  qt  K.  q 
ij  J  i 


(6 . 28) 


which  are  related  to  the  coefficients  c^ .  via  Eqs.  (6.21),  (6.18) 
and  (6.13)  : 


c..  *  (q.  K.  q)  a.  =  e..  a. 
ij  J  i  i  ij  i 


(6.29) 


Therefore  the  optimality  criterion  stated  in  Eq.  (6.27)  takes 
again  the  customary  "energetic"  form  : 


E  r  .  c  .  .  =  c' 
J  J  ij 


(6.30) 


where  the  e.j's  are  virtual  energy  densities  per  unit  weight  : 


78 


e 


i  j 


U_ 

.  a  . 
1  1 


(6.31) 


Note  that  the  design  variables  a^  are  still  defined  as  the 

transverse  sizes  of  the  structural  members,  so  that  P.l.a. 

th  ill 

represents  the  weight  of  the  i  element. 


In  the  special  case  where  only  one  flexibility  constraint  is 
specified,  the  optimality  criterion  states  that  the  virtual 
strain  energy  density  must  be  the  same  in  each  element.  The 
optimality  criterion  can  then  be  employed  to  generate  analy¬ 
tically  the  solution  to  the  explicit  problem  (6.24-6.26).  By 
following  out  the  same  calculations  as  in  section  2.2.2  of 
Ref.  [  1  ]  ,  the  optimum  design  variables  can  be  expressed  ex¬ 
plicitly  in  terms  of  known  quantities.  Taking  Eq.  (6.27)  with 

only  one  lagrangian  multipliers  (m=  I  )  yields 

1 

c  .  p+  1 

ai  -  (r  p  7-j-)  i  =  I  ,n  (6.32) 

'  i  i 


where  n  is  the  number  of  active  design  variables.  Substituting 
then  this  latter  equation  into  the  constraint  equation 

% 

n  c . 


i£. 


l 


u 


u 

o 


(6.33) 


l 


(where  is  the  contribution  of  the  passive  design  variables) 
leads  to  an  explicit  expression  for  the  lagrangian  multiplier  r 

P  ^  1  p 

1  J,  I",.  c.p] 


(r  p)ptl 


u  -  u 


[  \  *k  ck' 


P+1 


(6.34) 


Finally,  by  reintroducing  Eq .  (6.34)  into  Eq .  (6.32),  it  follows 
that 


a  . 
l 


{  — 
u 


n 


1  _J_ 

iP  ( _ LsP+  1 

‘  '  n .  I  .  ' 
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This  latter  equation  gives  the  solution  of  the  explicit  problem 
provided  the  subdivision  of  the  design  variables  into  active  and 
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passive  groups  is  correctly  defined.  It  is  worthwhile  mentioning 
that  Eq .  (6.35)  is  well  suited  for  the  design  of  plates  in  ben¬ 

ding.  Since  then  p  =  3,  the  redesign  relation  (6.35)  involves 
the  fourth  root  of  the  flexibility  coefficients  c.,  rather  than 
the  third  root  as  employed  by  ARMAND  and  LODIER  [  75  ]  on  an 
intui tive  basis  . 
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7.  COMPUTER  PROGRAM  IMPLEMENTATION 


This  chapter  gives  some  indications  about  the  computer  im¬ 
plementation  of  the  concepts  developed  in  the  present  work. 

The  resulting  optimization  capabilities  are  fully  integrated 
in  the  general  purpose  finite  element  program  SAMCEF  (15  ]  , 
which  is  applicable  to  large  structural  systems  (several  thou¬ 
sands  of  degrees  of  freedom  and  finite  elements).  The  optimi¬ 
zation  module  of  SAMCEF  is  built  to  loop  on  the  general  static, 
dynamic  and  stability  analysis  modules.  This  implies  that  all 
the  possibilites  offered  by  these  modules  are  still  available, 
as  well  as  those  of  the  auxiliary  modules  like  mesh  generators, 
plotting  modules  for  input  and  output,  etc...  Given  a  finite 
element  model,  the  user  may  ask  for  one  or  more  optimization 
steps,  without  any  thing  else  to  do  than  to  define  the  constraints. 

At  present  time  the  elements  whose  dimensions  are  taken  into 
account  in  the  operational  version  of  the  program  are  limited 
to  the  bars  and  plates  elements  in  extension.  The  next  version 
will  include  bending  elements  (beams  and  flat  shells)  according 
to  the  indications  given  in  chapter  6.  However  the  structure 
may  be  idealized  using  any  other  type  of  element  but  their  di¬ 
mensions  remain  unchanged  by  the  optimization  process.  In  par¬ 
ticular  the  possibility  of  using  super-elements  exists  and  is 
very  useful  for  representing  pre-op timized  or  fixed  parts  in  the 
structure.  In  absence  of  specifications  the  thickness  or  cross- 
section  of  each  finite  element  is  taken  as  a  design  variable  and 
is  allowed  to  be  i ndepe nden t ly  r e s i ze d .  However  the  finite  elements 
can  be  grouped  in  such  a  way  that  one  design  variable  is  assigned 
to  each  group.  The  possibility  exists  in  membrane  elements  to 
represent  composite  materials  like  reinforced  resins  as  the  su- 
•  rpiii tion  of  a  number  of  layers  with  independent  ortho  tropic 
rties.  The  thickness  of  each  layer  is  then  a  separate  de  - 
:  i  in  lc  so  that  the  superposition  of  results  allows  for 
••  ’  tion  of  the  composite. 
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In  the  current  version  of  the  SAMCEF  program,  three  distinct 
optimization  modules  are  provided  for  treating  separately  static 
constraints  on  stresses  and  displacements  (see  chapter  3),  dy¬ 
namic  constraints  on  natural  frequencies  (see  chapter  4 )  and 
stability  constraints  on  linear  buckling  loads  (see  chapter  5). 

It  is  envisioned  that  the  next  version  will  include  the  possibility 
of  taking  simultaneously  into  account  any  combinations  of  these 
three  types  of  behavior  constraints  by  using  a  data  base  system. 

The  program  can  evidently  make  use  of  the  restart  capabilities 
of  the  analysis  modules  but  it  is  also  designed  for  an  interactive 
use  when  the  designer  examines  the  solution  after  each  optimization 
stage.  After  each  structural  analysis  the  program  can  be  stopped, 
and  then  automatically  restarted  without  repeating  the  analysis. 

The  selection  of  the  optimization  algorithm  is  left  to  the 
user  but  recommendations  based  on  examples  are  provided  to  help 
the  choice  which  will  however  probably  remain  problem  dependent. 

7.1.  Fabrication  requirements 

In  practical  structural  design  problems  it  is  necessary 
to  take  into  account  various  f abricational  limitations.  In 
addition  to  the  usual  side  constraints,  which  allow  the  designer 
to  prescribe  realistic  minimum  or  maximum  gauge  sizes,  design 
variable  linking  is  often  necessary  to  permit  generation  of  a 
practically  meaningful  design.  Also  the  computer  program  must 
include  the  possibility  of  treating  fixed  members. 

In  order  to  insure  sufficient  accuracy  in  the  structural 
analysis,  finite  element  models  usually  involve  large  numbers 
of  elements  but  a  relatively  small  number  of  independent  design 
variables,  because  of  the  subdivision  of  the  various  members 
into  some  preselected  groups  (see  section  1.4).  As  expressed  in 
Kq .  (1.3)  design  variable  linking  consists  of  equality  constraints 

on  the  member  sizes  and  it  can  therefore  be  easily  handled  in  the 
problem  formulation. 
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Using  the  convenient  notation  [  74  ] 
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where  k  G  i  implies  "all  elements  that  are  linked  to  the  i 
independent  design  variable",  the  explicit  problem  statement 


takes  the  form 

(see 

Eqs.  2.18-2.20)  : 
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Note  that  l '  represents  the  cumulated  weight  of  all  elements 
in  the  group  when  a^  =  I.  When  some  or  all  stress  cons¬ 

traints  are  treated  using  the  FSD  procedure  (zero  order  appro¬ 
ximation,  see  section  3.2),  they  are  included  in  the  minimum 
size  limits  using  the  modified  stress  ratio  formula 
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From  the  point  of  view  of  the  minimization  algorithms,  no  change 
has  to  be  introduced  in  the  primal  and  dual  solution  schemes 
discussed  in  chapter  2,  because  the  explicit  problem  after  lin¬ 
king  exhibits  exactly  the  same  form  as  the  explicit  problem 


before  linking.  The  only  effect  of  linking  is  to  reduce  the 
number  of  independent  design  variables  from  n  to  n'. 

For  various  reasons  it  is  sometimes  interesting  to  incor¬ 
porate  in  the  structural  model,  finite  elements  that  are  not 
associated  with  any  design  variable  (superelements,  fixed  mem¬ 
bers,  etc...).  The  stiffness  and  mass  properties  of  such  ele¬ 
ments  are  not  affected  by  redesign,  however  their  contributions 
to  the  approximate  constraints  statement  must  be  taken  into  ac¬ 
count  in  the  form  of  fixed  terms  to  be  substracted  from  the  upper 
bounds  u..  For  example,  when  the  virtual  load  technique  is  em¬ 
ployed  to  generate  explicit  approximations  of  flexibility  cons¬ 
traints,  Eqs.  (3.7,  3.8)  should  be  rewritten 
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where  n'  is  the  number  of  variables  and  K.  are  the  stiffness 

1  o 

matrices  of  the  fixed  members  (i.e.  for  i  >  n'). 


7 . 2  Optimization  algorithms 

Four  optimization  algorithms  are  available  in  the  SAMCEF 
program.  The  user  can  select  any  one  of  them  depending  upon 
the  characteristics  of  each  specific  problem  :  the  number  o f 
independent  design  variables,  the  number  of  behavior  constraints 
and  the  expected  degree  of  non  linearity  of  the  constraints. 

7.2.1  PRIMAL  I  optimizer 

PRIMAL  1  is  a  first  order  projection  algorithm  based  on  the 
well  known  gradient  projection  method  for  linear  constraints. 
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It  uses  an  orthogonal  projection  operator  to  generate  a  sequence 
of  search  directions  that  are  constrained  to  reside  in  the  sub¬ 
space  defined  by  the  set  of  active  constraint  hyperplanes.  The 
successive  search  directions  are  conjugated  to  each  other  as 
long  as  there  is  no  change  in  the  set  of  active  constraint. 

The  PRIMAL  1  optimizer  operates  in  the  space  of  the  reciprocal 
design  variables  and  it  produces  a  sequence  of  steadily  improved 
feasible  designs  with  respect  to  the  linearized  problem.  Hence 
PRIMAL  I  can  be  adequately  used  for  seeking  a  partial  solution 
to  each  linearized  problem,  in  such  a  way  that  the  constraints 
of  the  primary  problem  remain  almost  satisfied.  This  is  achieved 
by  prescribing  an  upper  limit  on  the  number  of  one -d imens iona 1 
minimizations  performed  before  updating  the  approximate  problem 
statement  (parameter  k  ;  see  section  2.3).  PRIMAL  1  is  thus  the 
recommended  option  when  the  constraints  of  the  primary  problem 
are  highly  nonlinear  in  the  reciprocal  variables  (strong  struc¬ 
tural  redundancy).  The  algorithm  is  described  in  detail  in 
section  6.1  of  Ref.  [  1  ]  . 

7.2.2  PRIMAL  2  optimizer 

PRIMAL  2  is  a  second  order  projection  algorithm  especially 
well  suited  to  the  solution  of  problems  with  separable  objective 
function  and  linear  constraints.  It  uses  a  weighed  projection 
operator  to  generate  a  sequence  of  Newton's  search  directions  in 
the  subspace  formed  by  the  intersections  of  the  active  constraint 
hyperplanes.  PRIMAL  2  exhibits  the  same  features  as  PRIMAL  I, 
but  it  is  far  more  efficient.  Therefore  the  PRIMAL  2  optimizer 
is  well  suited  to  solve  exactly  each  linearized  problem,  in  which 
case  it  produces  the  same  iteration  history  as  the  dual  methods. 

PRIMAL  2  is  thus  a  recommended  option  when  the  behavior  cons¬ 
traints  are  very  shallow  in  the  space  of  the  reciprocal  variables 
(weak  structural  redundancy).  Note  however  that  the  DUAL  2  option 
is  usually  more  efficient.  The  PRIMAL  2  algorithm  is  described 
in  section  6.2  of  Ref.  (I  ]. 
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7.2.3  DUAL  2  optimizer 

The  dual  method  formulation,  which  exploits  the  separable 
form  of  the  approximate  problem,  consists  in  maximizing  the 
explicit  dual  function  subject  to  nonnegativity  constraints 
on  the  dual  variables.  As  explained  in  section  2.4.1,  this 
approach  is  very  efficient,  because  the  dimensionality  of  the 
dual  space  is  primarily  dependent  on  the  number  of  critical 
behavior  constraints,  which  is  relatively  low  for  many  struc¬ 
tural  optimization  problems  of  practical  interest. 

DUAL  2  is  a  dual  method  which  employs  a  second  order  Newton 
type  of  algorithm  to  find  the  maximum  of  the  dual  function  when 
all  the  design  variables  are  continuous.  It  operates  in  a  se¬ 
quence  of  dual  subspaces wi th  gradually  increasing  dimensions, 
so  that  the  effective  dimensionality  of  the  dual  problem  does 
not  exceed  the  number  of  active  behavior  constraints  by  more 
than  one.  Since  the  DUAL  2  optimizer  has  been  found  to  be  highly 
efficient  in  practice,  it  is  the  recommended  option  for  pure 
continuous  variable  problems,  unless  the  number  of  active  behavior 
constraints  is  expected  to  be  high.  The  algorithm  is  described 
in  section  4.3  of  Ref.  [1  ]  ,  as  well  as  in  section  3  of  Ref.  |  I  3]  . 

7.2.4  DUAL _ 1  optimizer 

DUAL  I  is  a  dual  method  which  employs  a  specially  devised  first 
order  gradient  projection  type  of  algorithm  to  find  the  maximum 
of  the  dual  function  when  all  or  some  of  the  design  variables  are 
discrete  (see  Ref.  [19  ]).  The  DUAL  1  algorithm  incorporates 
special  features  for  handling  the  dual  function  gradient  discon¬ 
tinuities  that  arise  from  the  primal  discrete  variables.  These 
discontinuities  occur  on  specific  hyperplanes  in  the  dual  space. 

The  DUAL  I  algorithm  determines  usable  search  directions  by  pre¬ 
lecting  the  dual  function  gradient  on  the  intersection  of  the 
successively  encountered  first  order  discontinuity  planes.  It 
should  be  noted  that  the  DUAL  I  optimizer  remains  applicable  to 
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pure  continuous  variables  problems,  in  which  case  it  reduces  to 
a  special  form  of  the  conjugate  gradient  method.  However  it  is 
generally  less  efficient  than  the  DUAL  2  optimizer.  DUAL  1  was 
initially  conceived  for  the  ACCESS  3  program  [  14  ]  and  it  is 
described  in  detail  in  chapter  4  of  Ref.  [  13  ]  . 

7.2.5  Core  requirements 

Because  of  their  special  implementation,  which  takes  advan¬ 
tage  from  the  simple  form  of  the  side  constraints,  the  projection 
algorithms  PRIMAL  1  and  PRIMAL  2  require  a  modest  core  size  (see 
Ref.  [1  ]).  The  number  of  words  necessary  to  solve  a  given  pro¬ 
blem,  with  n  independent  design  variables  and  n  linearized  be¬ 
havior  constraints,  is  given  by  • 

nbr  of  words  =  n  x  (m  +  10)  +  — — — — —  -  +  3  x  m  (7.10) 

This  formula  also  applies  to  the  DUAL  2  optimizer,  however  it 
should  be  modified  for  the  DUAL  1  optimizer  when  discrete  va- 
variables  are  in vo 1 ve d . F i g  .  5  represents  graphically  the  core  re¬ 
quirement  given  by  Eq.  (7.10).  It  should  be  clearly  recognized 
that,  because  design  variable  linking  is  most  often  employed 
in  practical  applications,  the  SAMCEF  program  is  capable  of 
dealing  with  structural  optimization  problems  involving  thousands 
of  finite  elements.  For  example  a  problem  with  500  design  va¬ 
riables  -  which  might  well  correspond  to  5000elements  -  and  20 
linearized  behavior  constraints  can  fit  in  a  core  of  less  than 
16000  words  (on  IBM  370-158  the  computer  program  requires  then 
256  K).  Note  that  no  core  limitation  is  associated  with  the 
analysis  modules  of  SAMCEF,  because  they  are  organized  in  such 
a  way  that  they  can  solve  very  large  problems  using  a  modest 
size  central  core.  They  employ  a  frontal  equation  solver  with 
substructing  and  extensive  peripheric  storage. 

It  can  be  seen  from  examining  Fig.  5  that  the  main  limitation 
arises  from  the  number  m  of  behavior  constraints  retained  in  the 


linearized  problem  statement.  If  the  number  of  constraints  is 
raised  up  to  100  in  the  previous  example  (with  500  design  vari¬ 
ables),  then  the  central  core  requirement  increases  to  60000 
words.  In  addition  it  is  evident  that  linearizing  the  behavior 
constraints  demands  a  large  computational  effort,  because  this 
implies  treating  additional  loading  cases  in  the  structural  re¬ 
analyses  (see  section  3.1).  Finally  the  computer  time  expanded 
in  the  optimizer  itself  can  become  prohibitive  when  the  number 
of  linearized  constraints  is  large.  This  is  apparent  in  the  case 
of  dual  algorithms,  since  the  dimensionality  of  the  dual  problem 
is  precisely  equal  to  the  number  of  linearized  constraints.  It 
is  also  true  in  the  case  of  primal  algorithms,  because  most  of 
the  computational  effort  quickly  increases wi th  the  number  of 
linear  constraints  (construction  of  the  projection  matrices, 
evaluation  of  the  maximum  allowable  step  length,  selection  of 
the  set  of  active  constraints,  etc...). 

Therefore  it  is  important  to  reduce  as  much  as  possible  the 
number  of  behavior  constraints  retained  at  each  stage  of  the 
optimization  process.  Constraint  diletion  techniques,  such  as 
those  proposed  by  SCHMIT  and  MIURA  [  6-8  ] ,  are  simple  but  ef¬ 
fective  means  of  achieving  that  goal.  In  this  approach  only  the 
critical  and  potentially  critical  behavior  constraint  are  included 
in  the  linearized  problem  statement  at  each  redesign  stage.  Po¬ 
tentially  critical  constraints  are  defined  as  those  which  are 
close  to  their  allowable  limit  within  a  given  tolerance  (which 
can  be  dynamically  updated  at  each  stage).  Note  that  the  struc¬ 
tural  reanalysis  must  be  decomposed  in  two  parts  :  first  the 
constraint  values  are  computed  for  the  real  loading  conditions, 
and  then  the  gradients  of  the  retained  constraints  are  evaluated 
using  virtual  or  pseudo-load  cases  (see  section  3.1).  Finally, 
using  zero  order  approximation  for  the  stress  constraints  also 
permits  a  dramatic  reduction  in  the  number  of  linearized  cons¬ 
traints  (see  section  3.2). 
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7 . 3  Program  organization 

On  the  basis  of  the  foregoing  developments,  we  propose  in 
this  section  the  program  organization  that  seems  the  most  appro¬ 
priate  for  dealing  with  large  structural  optimization  problems  : 

( I )  Preprocessor 

-  generation  of  the  element  stiffness  matrices,  mass 
matrices  and  implicit  loads  for  unit  values  of  the 
design  variables 

-  data  preparation  depending  on  design  variable  lin¬ 
king 

loop  2 )  Update  of  the  element  stiffness  matrices,  mass  matrices 
and  implicit  loads  for  the  current  values  of  the  design 
variables 

(  3  )  Static  ana ly s i s 

-  triangular  decomposition  of  the  system  stiffness 
matrix 

-  evaluation  of  the  displacements  (back  substitutions 
for  the  real  loading  cases) 

-  evaluation  of  the  stresses 

( 4 )  Dynamic  analysis 


-  assemblage  of  the  system  mass  and  stiffness  matrices 

-  solution  of  the  e i genprob 1 em ,  yielding  the  frequencies 
and  the  associated  modal  displacements 

(  3 )  Stability  analysis 

-  generation  of  the  element  geometric  stiffness  matrices 
(from  the  stresses  previously  computed) 


i 
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-  assemblage  of  the  system  matrices  and  solution 
of  the  eigenproblem 


(6)  Seal ing  of  the  design  variables  and  convergence 
checks  - - >  (10) 


(  7 )  Linearized  problem  statement 

-  selection  of  the  potentially  critical  constraints 

-  selection  of  zero/first  order  approximation  scheme 
for  each  stress  constraint 

-  definition  of  the  additional  dummy  loads  (virtual 
or  pseudo  load  cases) 

-  evaluation  of  the  constraint  gradients  (back  sub¬ 
stitutions  for  the  dummy  load  cases) 

( 8 )  Optimi zation 

-  solution  of  the  linearized  problem  using  either 
primal  or  dual  algorithm 

1  oop ( 9 )  S  torage  of  the  design  variable  vector  and  return  to  (2) 


(10)  Print  out  of  the  results 


End  . 


In  this  computer  program  organization,  it  is  assumed  that  the 
search  for  the  optimal  design  is  performed  in  a  fully  automatic 
way.  It  should  be  noted,  however,  that  it  is  often  desirable  to 
allow  for  a  human  intervention  in  the  redesign  process.  For 
example  very  large  structural  models,  because  they  are  time  con¬ 
suming,  cannot  be  treated  in  a  single  run  and  it  is  better  exe¬ 
cuting  the  optimization  program  stage  by  stage  with  intermediate 
verification  of  the  results.  The  control  parameters  can  then  be 
reset  periodically  to  adequate  values  (change  of  the  optimization 
algorithm,  modification  of  the  £  parameter,  specification  of  new 
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tolerances  for  the  selection  of  the  linearized  constraints, 
etc...).  By  storing  in  a  data  base  all  the  analysis  results 
at  each  stage,  it  is  moreover  possible  to  restart  the  program 
at  any  previously  generated  design  point. 

In  this  connection  ,  the  concept  of  automatic  redesign 
should  be  replaced  by  that  of  interactive  redesign  on  a  gra¬ 
phic  terminal,  allowing  therefore  the  designer  to  easily  mo¬ 
nitor  the  optimization  process.  It  is  worth  pointing  out  that 
the  high  quality  explicit  approximations  of  the  behavior  cons¬ 
traints  employed  in  this  work  could  also  be  used  in  the  context 
of  computer  aided  design  procedures.  Indeed  the  constraint 
gradients  with  respect  to  the  reciprocal  design  variables  can 
be  viewed  as  "sensibility  indices"  whose  values  provide  quali¬ 
tative  and  even  quantitative  informations  about  the  way  in 
which  the  structure  should  be  redesigned. 


8.  NUMERICAL  EXAMPLES 


A  wide  variety  of  truss  and  box  structures  were  proposed  in 
Ref.  [  1  ]  to  illustrate  the  power  of  the  mixed  method  and  gene¬ 
ralized  optimality  criterion  concepts  applied  to  stress  and  dis¬ 
placement  constraints.  In  this  section  attention  is  focused  on 
sample  problems  involving  frequency  and  buckling  constraints,  as 
well  as  problems  with  beam  elements. 

Although  the  results  presented  in  this  chapter  are  only  con¬ 
cerned  with  the  SAMCEF  program  [  15  ]  ,  it  should  be  emphasized  that 
other  significant  problems  have  been  successfully  treated  with  the 
ACCESS-3  program  [  13,  14  ]  ,  in  which  the  dual  method  approach  is 

also  available. 

8.1.  E  las t i c  r od  w i th  tip  mass 


In  order  to  validate  the  optimization  strategy  proposed  in  sec¬ 
tion  4.3,  the  following  simple  problem  was  considered  :  minimize 
the  weight  of  an  elastic  bar  subject  to  minimum  fundament al  frequen¬ 
cy  constraint  (longitudinal  vibration).  The  bar  is  fixed  at  one 
end  and  it  supports  a  given  mass  at  the  free  end  (see  Fig.  6)  . 

This  example  is  interesting  because  it  has  received  an  exact 
analytical  solution  in  the  context  of  continuous  models  studies 
[  34,  35,  38  ] .  The  optimal  structure  for  a  prescribed  lower  limit 
a)  is  defined  by  the  cross-sectional  area  a(x)  and  the  axial  displa¬ 
cement  u(x)  as  follows  : 

a  ( x)  =  th  (tU)  — -  (8.1) 

P  cti  (  8  x) 


u  (x) 


sh ( 6x) 
shtlU 


(8.2) 
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In  these  expressions,  />  denotes  the  mass  density  for  the  given 
material,  E  is  its  Young  modulus,  l  represents  the  length  of  the 
rod  and  m  is  the  tip  mass  (at  x  =  £).  The  total  mass  for  the  op¬ 
timal  structure  (m  non  included)  is 

W*  =  in  sh2(B£)  (8.4) 

These  analytical  results  were  compared  with  those  generated  by 
using  recursively  the  optimality  criterion  stated  in  Eq .  (4.36). 

For  a  finite  element  model  made  up  of  10  bars,  an  optimal  weight 
of  0.6647  kg  is  produced  in  8  iterations.  The  corresponding  ana¬ 
lytical  solution  (8.4)  leads  to  W*  =  0.6646  kg.  The  numerical 
values  of  the  cross-sectional  areas  and  the  modal  displacements 
in  the  finite  element  solution  are  almost  identical  to  the  exact 
values  provided  by  Eqs.  (8.1)  and  (8.2)  (up  to  the  fourth  digit). The 
iteration  history  for  this  10-bar  model  is  illustrated  in  Fig.  7. 

The  same  problem  was  also  treated  with  a  cruder  finite  element 
model  involving  only  4  bars.  The  optimum  is  attained  in  9  itera¬ 
tions.  Table  I  shows  the  detailed  iteration  history.  The  final 
design  weighs  0.6656  kg  and  it  remains  close  to  the  analytical 
solution  despite  the  simplicity  of  the  discretization  (only  4  de¬ 
grees  of  freedom).  The  variations  of  the  cross-sectional  area  and 
of  the  axial  displacement  are  depicted  in  Fig.  6.  The  agreement 
between  the  numerical  results  provided  by  the  optimality  criterion 
(4.36)  and  the  analytical  solution  given  by  Eqs.  (8.1,  8.2)  is  ex¬ 
cellent. 

8.2.  Sandwich  beam  with  Euler  buckling  contraint 

Attention  is  now  directed  to  the  Euler  sandwich  column  shown 
in  Fig.  8.  The  beam  has  length  l  and  is  hinged  at  both  ends.  The 
problem  consists  in  minimizing  its  weight  for  a  given  critical  bu¬ 
ckling  load  P  (see  Fig.  8).  For  a  sandwich  beam  with  rectangular 
cross-section,  specific  weight  and  stiffness  properties  (per  unit 
length)  are  linearly  related,  so  that  the  problem  can  be  treated 
by  using  the  same  strategies  as  in  the  case  of  thin-walled  structures. 
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Considering  thus  a  rectangular  sandwich  beam  with  constant 
core  height  2h  and  constant  width  b  and  assuming  identical  upper 
and  lower  face  sheets  with  variable  thickness  t  «  h,  it  follows 
that  the  specific  bending  stiffness  can  be  written 


s  =  2  E  b  h  t 


(8.5) 


and  that  the  flange  specific  weight  takes  the  form 


w  =  2  p  b  t 


(8.6) 


where  E  and  p  represent  respectively  the  Young's  modulus  and  the 
mass  density  for  the  flange  material.  Note  that  only  the  contri¬ 
bution  of  the  face  sheets  to  the  structural  weight  can  be  minimized, 
because  the  core  sizes  are  kept  constant.  The  onjective  function 
is  therefore  linear  in  the  design  variable  t  (sheet  thickness). 


The  analytical  solution  to  this  problem  was  obtained  by  PRAGKR 
and  TAYLOR  [  38  ] .  At  the  optimum  the  specific  stiffness  s(x)  and 
the  lateral  displacement  u(x)  that  characterizes  the  critical  bu¬ 
ckling  mode  are  expressed  as  follows  : 


s(x)  =  |  x  («.  -  x)  (8.7) 

u  ( x )  =  x  (Z  -  x)  (8.8) 

a 


The  sheet  thickness  can  then  be  evaluated  using  Eq  .  (8.5)  and 
Eq .  (8.7).  Integrating  the  specific  weight  stated  in  Eq .  (8.6) 
along  the  length  of  the  column  yields  the  optimal  weight 


W 


x 


P  p  9  3 
12  E  h2 


(8.9) 


If  W°  denotes  the  weight  of  a  beam  with  uniform  thickness,  sup¬ 
porting  the  same  buckling  load  P,  it  follows  that 


W 


M 


W°  =  0.822  W 
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The  latter  expression  demonstrates  the  benefit  gained  from  opti¬ 
mizing  the  column. 


Using  a  finite  element  model  involving  10  equal  length  segments, 
the  discretized  solution  resulting  from  applying  the  optimality 
criterion  redesign  relations  (see  Eq .  5.29)  is  generated  in  6  ite¬ 
rations.  Table  2  reproduces  the  iteration  history.  Note  that  by 
symmetry  only  5  beam  elements  are  necessary  to  describe  the  problem. 
The  final  design  weighs  265.7  kg,  while  the  analytical  solution  (8.9) 
leads  to  260  kg.  The  continuous  and  discretized  solutions  are  com¬ 
pared  in  Fig.  8  (see  also  Table  2) .  It  can  be  seen  that  the  nume¬ 
rical  results  are  very  close  to  the  analytical  solution,  although 
the  finite  element  model  involves  only  10  degrees  of  freedom. 


8.3.  Euler  column  with  rectangular  cross-section 

It  has  been  shown  in  section  6.1.1  that  the  optimum  design  of 
beam  structures  depends  upon  the  relation  between  the  moment  of 
inertia  and  the  cross-sectional  area  of  the  beam  members.  In  order 
to  illustrate  the  various  possibilities  permitted  by  Eq .  (6.1),  we 

consider  again  the  Euler  column  problem,  but  we  assume  now  a  rec¬ 
tangular  cross-section  with  height  h  and  width  b.  As  explained  in 
section  6.1.1,  the  case  where  b  is  variable  and  h  is  constant  reduces 
to  the  previously  examined  problem,  because  the  bending  stiffness 
for  each  member  is  linear  with  respect  to  the  cross-sectional  area. 
Now,  if  b  is  kept  constant  and  h  is  variable,  the  relation  (6.5) 
holds,  in  which  case  one  should  choose  p  =  3  in  the  expression  (6.8) 
of  the  stiffness  matrix.  Finally  if  the  cross-section  shape  remains 
constant  along  the  beam  axis  (i.e.  fixed  h/b  ratio)  one  should  adopt 
p  =  2  as  shown  by  Eq .  (6.6). 


From  the  following  remarks  it  follows  that  the  optimum  design 
must  be  sought  by  applying  the  optimality  criterion  redesign  rela¬ 
tions  (6.27)  with  p  =  1,  2  or  3  depending  upon  the  hypothesis  made 
on  the  variable  quantity.  Of  course  it  should  be  kept  in  mind  that 
u  represents  here  the  critical  buckling  load  instead  of  a  flexibi¬ 


lity 


Consequently  the  coefficients  c^j  must  be  corrected  according 


to  Eq  .  (5.25),  in  order  to  take  the  scaling  effect  into  conside¬ 

ration.  Since  only  one  behavior  constraint  is  imposed  in  the 
Euler  column  problem  statement,  the  redesign  relation  (6.35)  is 
directly  applicable  once  the  adequate  value  of  p  is  selected. 

In  each  of  the  three  cases  considered  the  finite  element 
model  involves  5  beam  members  (symmetry).  In  the  initial  design 
the  width  and  the  height  of  the  rectangular  cross-section  are 
taken  as  10  and  I,  respectively,  which  leads  to  a  weight  equal  to 
10  (arbitrary  unit  system  ;  see  Fig.  9).  After  scaling  of  the 
design  variables,  the  weight  of  the  strictly  critical,  feasible 
design  takes  on  different  values  in  the  three  cases  considered, 
because  the  scaling  factor  multiplies  the  moments  of  inertia  of 
each  beam  element,  rather  than  their  cross-sectional  areas. 

Table  3. a  contains  the  iteration  histories  for  each  three  si¬ 
tuations  p  =  1,  p  =  2  and  p  =  3.  Only  four  to  six  structural  re¬ 
analyses  are  required  to  achieve  convergence.  Table  3.b  represents 
the  corresponding  final  designs,  which  all  satisfy  the  optimality 
criterion  stated  in  Eq .  (6.30).  The  least  weight  design  is  obtained 
when  the  width  of  the  cross-section  is  varied  (p  =  I).  From  the 
tip  to  the  middle  of  the  beam,  the  moment  of  inertia  (and  thus  the 
cross-sectional  area)  is  increased  by  a  factor  of  five.  The  hea¬ 
viest  design  is  obtained  by  choosing  the  height  of  the  cross-sec¬ 
tion  as  the  variable  quantity  (p  =  3).  Now  the  moment  of  inertia 
is  modified  by  a  factor  of  ten  along  the  beam.  However,  the  change 
in  cross-sectional  area  is  only  by  a  factor  of  two,  which  explains 
the  higher  value  of  the  optimal  weight  in  this  case.  The  variations 
of  tiie  moment  of  inertia  and  of  the  cross-sectional  area  are  dis¬ 
played  in  Fig.  9  for  the  three  cases  considered. 


The  sample  problems  examined  in  the  previous  sections  are  con¬ 
cerned  with  a  single  behavior  constraint,  in  which  case  each  expli¬ 
cit  problem  generated  in  sequence  can  be  solved  in  closed  form 
(see  for  example  section  6.3).  When  multiple  behavior  constraints 
are  simultaneously  active,  application  of  the  generalized  optimality 
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criterion  stated  in  Eqs.  (2.21,  2.25)  requires  the  explicit  pro¬ 
blems  to  be  solved  by  using  dual  algorithms.  As  an  example  we 
consider  the  I-beam  represented  in  Fig.  10.  Design  data  are  gi¬ 
ven  in  Table  4,  including  the  lower  and  upper  limits  imposed  on 
the  natural  frequencies. 

The  problem  consists  in  minimizing  the  weight  of  the  beam 
while  controlling  the  frequencies  of  its  three  first  eigenmodes  : 

.  flange  flexion  (along  the  Y-axis) 

.  torsion  (with  respect  to  the  X-axis) 

.  web  flexion  (along  the  Z-axis) 

From  the  structural  analysis  point  of  view,  a  difficulty  in  this 
problem  is  to  find  a  model  made  up  of  membrane  elements,  that  is 
capable  of  representing  properly  each  three  vibration  modes.  The 
two  bending  modes  are  well  taken  into  account  in  a  membrane  model. 
However  the  torsional  mode  is  not  correctly  represented,  because 
it  reduces  to  differential  flexion  of  the  two  flanges  (with  the 
same  frequency  as  the  fundamental  mode).  It  is  therefore  necessary 
to  include  fictitious  diaphragms  in  the  analysis  model,  so  as  to 
prevent  differential  flange  flexion.  These  diaphragms  are  intro¬ 
duced  with  zero  mass  density. 

Table  5  reproduces  the  natural  frequencies  of  the  three  first  eigen¬ 
modes,  evaluated  on  one  hand, with  pure  membrane  models  (with  and 
without  diaphragms)  and,  on  the  other  hand,  with  a  much  more  accu¬ 
rate  model  made  up  of  flat  shell  elements.  In  both  cases,  the 
beam  is  subdivided  into  five  equal  length  segments,  which  leads  to 
25  rectangular  elements  (plus  10  diaphragms  in  the  third  model  con¬ 
sidered  in  Table  5).  The  model  retained  yields  excellent  results, 
which  are  in  good  agreement  with  those  generated  using  flat  shell 
elements.  It  involves  250  degrees  of  freedom,  because  a  second 
degree  displacement  field  is  taken  in  each  membrane  element.  In 
the  numerical  solution  of  the  eigenproblem  all  the  interface  degrees 
of  freedom  are  condensed,  so  that  only  87  degrees  of  freedom  are 
retained.  It  should  be  noted  that  in  the  optimization  process  the 
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thickness  of  the  diaphragms  are  fixed  to  their  initial  value 
(I  mm),  so  that  the  problem  involves  25  design  variables. 

8.4.1.  Case  A  ;  inequality  constraints  on  the  natural  fre¬ 
quencies 

In  a  first  optimization  exercise,  lower  and  upper  bounds  are 
imposed  on  each  three  eigenf requencies  (see  Table  3).  Table  6 
and  Fig.  II  summarize  the  iteration  histories  by  reproducing  the 
variations  of  the  weight  and  frequencies  with  the  number  of  struc¬ 
tural  reanalyses.  The  convergence  of  the  optimization  process  is 
remarkably  rapid.  Only  5  analyses  are  sufficient  to  generate  an 
optimum  design,  weighing  393.7  kg.  The  fundamental  frequency  does 
not  reach  the  prescribed  lower  bound  (1.0  Hz)  nor  the  upper  bound 
(1.2  Hz).  The  two  higher  frequencies  are  equal  to  their  respective 
minimal  allowable  value  (1.2  Hz  and  2.5  Hz). 

The  final  design  is  given  in  Table  7.  Note  chat  the  thickness 
in  the  upper  flange  differs  from  those  in  the  lower  flange,  because 
of  the  support  condition  at  the  fifth  of  the  beam  length.  Also  an 
increase  in  the  flanges  thicknesses  appears  at  the  support  level. 
With  regard  to  the  web  thickness  it  can  be  seen  that  it  is  lower 
than  the  flange  thickness,  mainly  at  the  hinged  end  neighborhood. 

8.4.2.  Case  B  ;  equality  constraints  on  the  natural  frequencies 

In  a  second  optimization  exercise,  equality  constraints  were 
assigned  to  each  three  frequencies,  in  order  to  reduce  the  funda¬ 
mental  frequency  to  1.0  Hz,  while  keeping  the  two  other  frequencies 
at  1.2  Hz  and  2.5  Hz,  respectively.  The  iteration  history  presented 
in  Table  6  and  illustrated  in  Fig.  II  again  demonstrates  the  remar¬ 
kable  efficiency  of  the  generalized  optimality  criterion.  The  in¬ 
crease  in  weight  after  the  first  r  design  stage  is  due  to  the  equa¬ 
lity  constraint  which  is  now  imposed  on  the  first  frequency.  Note 
that  only  6  structural  reanalyses  are  required  to  achieve  conver¬ 
gence  and  that  the  equality  constraints  are  satisfied  within  less 
-4 

than  10  relative  accuracy. 
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Table  7  represents  the  final  designs.  Although  the  weights 
btained  in  the  two  cases  considered  are  very  close  to  each  other 
(393.7  kg  in  case  A  and  394.1  kg  in  case  B) ,  the  material  distri¬ 
bution  is  quite  different  in  case  B,  because  of  the  decrease  in 
the  fundamental  frequency.  The  differences  in  thickness  between 
the  two  flanges  are  more  pronounced.  Moreover  the  web  thickness 
varies  more  strongly  along  the  beam  axis. 


-J 
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EULER  COLUMN  WITH  RECTANGULAR  CROSS-SECTION 


width 


TABLE  4 

DESIGN  DATA  FOR  I-BEAM  PROBLEM 


Material 
Young's  modulus 
Poisson's  ratio 
mass  density 
minimum  thickness 
initial  thickness 
non-s tructural  mass 


steel 

E  =  2 . I05  N/mm2 
v  =  0.3 

p  *  7.8  10  ^  kg/mm^ 

£  =  1  mm 
a  =  10  mm 

20.000  kg  (see  Fig.  10) 


■ 

mode 

frequency 

limits  (Hz) 

type  of 

cons  traint 

minima  1 

maximal 

■1 

1 

1  .0 

1  .  2 

2 

1  .2 

2.5 

inequal i ty 

■ 

3 

2.5 

/ 

1 

I  .0 

1  .0 

Case  B 

2 

1  .2 

1  .  2 

equality 

3 

2.5 

2.5 

-  112- 

TABLE  5 

COMPARISON  OF  VARIOUS  ANALYSIS  MODELS  FOR  I-BEAM 


mode 

membrane 

membrane  2nd  degree 

flat  shell 

1st  degree 

without  diaphragm 

with  diaphragm 

2nd  degree 

1  flexion  (Y) 

l  .9782 

0.9285 

0.9285 

0.9286 

2  torsion  (X) 

1 .9782 

0.9285 

1  .0512 

1 .0379 

3  flexion  (Z) 

2.9515 

2.4970 

2.4970 

2.4951 

TABLE  6 

ITERATION  HISTORY  FOR  I-BEAM  PROBLEM 


ase  A  :  inequality  constraints!  case  B  :  equality  constraint! 


Iteration  we 
( 


frequencies  (Hz) 


507.00  0.9285  1.051  2.497 
443.14  1.137  1.241  2.614 
397.99  1.102  1.202  2.514 
394.11  1.101  1.200  2.501 
393.74  I. 103  1.200  2.500 


weight  frequencies  (Hz) 
(kg)  1 


507.00 

565.56 

432.78 

397.17 

394.06 

393.69 


0.9285 
1  .003 
0.9849 
1  .003 
1 .000 
1 .000 


2 

3 

1 

1.051 

2.497 

1.253 

2.679 

1.250 

2.531 

1.205 

2.508 

1.200 

2.501 

1.200 

2.500 

free 


THE  LINEARIZED  PROBLEM  CONVERGENCE  OF  THE  MIXED  METHOD 

FIG.  1  FIG.  2 


1 2L  CONVENTIONAL  kdi7  GENERALIZED 

_0  OPTIMALITY  1  . ,  r.  a0  OPTIMALITY 

U=U  oP  CRITERIA  \  7  CRITERION 


COMPARISON  OF  FIRST  AND  ZERO  ORDER  APPROXIMATIONS 


CONTINUOUS  DISCRETIZED 

MODEL  MODEL 


Vx 


young’s  Mooutus  E=  2  10”  N/m2 
MASS  DENSITY  p=  7800  kg/m3 
LENGTH  1=  Irn 

tip  mass  m  =  1kg  4 

INITIAL  CROSS-SECTIONAL.  AREA  a':  B  m2 
MINIMUM  FREQUENCY  ¥=  600  Hz 


FIG.  6  ELASTIC  ROD  WITH  TIP  MASS 


FIG.  7  ITERATION  HISTORY  FOR  ELASTIC  ROD 


BUCKLING  CONSTRAINT  EULER  COLUMN  WITH  RECTANGULAR  CROSS-SECTION 


